用EDA处理蛋白质数据小记

这两天在做蛋白质数据的EDA（exploratory data analysis）试了Table One、PCA、Volcano plot和Manhattan plot。实现语言为R。

数据简要介绍
蛋白质在全血中浓度log后的值，以及一张含各variable的表。protein数据648（proteins）*92（patients）。数据lung cancer case:control=322:325基本持平。case指lung cancer。

Table One
Table One用lung cancer case作为输出求各variable的p-value。Table One只能提供数据大致的第一印象。

需要注意的是在Table One中显著的variable并不能作为后面denoise data所用的adjust variable，原因有：

1. 有的variable间本来可能就有联系，比如BMI和高血压，如果都adjust可能会有重复调整的风险
2. Linear Mixed Model选择adjust variable应该根据reference选择合理的variable，是为了观察除这两个variable外潜在相关的variable，而不是根据p-value选择。

PCA
PCA如下图，结果很不好。问了别人，别人的结果也无法清晰分类。最后，大家的结论是PCA无法很好解释该数据，选择不用PCA。

Volcano Plot
因为笔者要探索smoking以外variable对lung cancer的影响，所以绘制了两张，一张before adjust by smoking, 一张after。每个点代表一个protein。纵坐标为-log10（adjust p-value），点越靠近图顶部表示差异越显著。横坐标为log2(fold change)，点越偏离中心，表示差异倍数越大。

名词解释

• fold change: 差异倍数，及数据在一组样品中的表达值的均值除以其在另一组样品中的表达值的均值。所以火山图只适合展示两组样品间的比较。
• 为什么做log2转化： fold change如果直接表示，值上到无限大，下至0，绘图很不方便。如果用log2，上调的数值都大于等于1，下调的数值都小于等于-1，绘图更方便。
  不过很遗憾，显著的点不多。
• 为什么adjust p-value: 统计学家认为多次的检测会引入假阳性结果，需要一个多重假设检验矫正。Bonferroni correction：每一次统计检验获得的p-value都乘以总的统计检验次数获得adjust p-value。 由于Bonferroni
  correction太过严厉，容易找不到有差异的数据，统计学家后来又提出别的相对宽松的计算方法，比如BH。BH获得的是假阳性率（FDR）。FDR<0.1表示我们对假阳性率的容忍度至多是10%。本次数值根据BH方法调整，有颜色代表significant。
• 为什么做-log10转换： 因为FDR是0-1/数值越小表示统计显著。-log10（adjust p-value）转换后正好是反过来，数值越大越显著。

Manhattan plot
将GWAS分析之后所有SNP位点的p-value在整个基因组上从左到右一次画出来，并且为了更直观的表达结果，通常会将p-value转换成-log10（p-value）。所以Y轴高度越高，关联性越强。
目前为止网上找到的都是用来画基因数据的。试着画了一下蛋白质数据，效果不好。

总之，univariate analysis并没找出什么有效protein。做完这些后会将数据用LMM denoise（adjust by age, gender, bmi）。

最后推荐一个无意间发现的宝藏网站http://www.ehbio.com/ImageGP/
很优秀的生信绘图网站。不过数据输入时格式要求有点多，有点麻烦。真的佩服这种交互功能。

Reference：

1. Volcano plot | 别再问我这为什么是火山图http://blog.sciencenet.cn/home.php?mod=space&uid=118204&do=blog&id=1222648
2. 如何理解GWAS中Manhattan plot和QQ plot所传递的信息。https://www.jianshu.com/p/987859ae503c