找差異蛋白
rm(list=ls()) #清空環(huán)境
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的數(shù)目惶楼，取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取對數(shù)
pvalue=NULL
#循環(huán) 將每個蛋白的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)
