找差异蛋白
rm(list=ls()) #清空环境
pcut=0.05
fccut=1.5
MaxNA=1 #生物学重复最大允许几个缺失
#第一列是蛋白名 2、3、4列是生物学重复
#把EXCEL保存为txt文件(另存为-制表符分割)
oridata=read.table("data/demo.pro.foldChange.txt",header=T,row.name=1)
#将0和10的蛋白表达量均替换为NA
oridata[oridata==0]=NA
oridata[oridata==10]=NA
#按行(蛋白)计算NA的数目,取0或者1个NA的蛋白
na_count=-apply(oridata,1,function(x)sum(is.na(x)))
matrix=oridata[na_count<=MaxNA,]
matrix=as.data.frame(matrix)
trData=log2(matrix) #对FC取对数
pvalue=NULL
#循环 将每个蛋白的log2(FC)与均值为0(即FC=1)做t检验,取p.value
for(i in (1:nrow(trData))){pvalue[i]=t.test(trData[i,],mu=0)$p.value}
FCmean=rowMeans(matrix,na.rm = T) #计算每个蛋白的FC均值
out=cbind(matrix,FCmean,pvalue) #输出结果包含原始矩阵、FC值、t检验的p值
write.csv(out,file="demo/out.all.csv",quote = F)
#abs是绝对值 log2X x大于1则是正值 X小于1则是负值
FCindex=abs(log2(FCmean))>abs(log2(fccut)) #取FC均值大于阈值的蛋白
pindex=pvalue<pcut
difindex=FCindex & pindex #dif表示对于FC和P值同时满足的
diff=out(difindex,) #选出满足的,付给新变量diff
write.csv(diff,file="demo.out.diff.csv",quote=F)