前言
“皇上,您还记得大明湖畔的夏雨荷吗?”
“小白,您还记得fam文件的第六列吗?”
当一个数据集(如肺癌),包括很多中类型的肺癌的时候
如腺癌(表型=2),鳞癌(表型=3),小细胞(表型=4)的时候,
你就需要把这部分数据提取出来了~
小白,不要告诉我你想在成千上万的数据中一个一个筛选出来~
首先,要达到“提取”=保留的目的,我们需要看需要什么指令:
这就是keep指令:
人话:相信我,小白这样看,一辈子都提取不出来
因为这里涉及到“准备中间的那个sampleID”文件。
所以,我们需要分两步走:
Step 1:弄这个“sampleID”文件
”sampleID”文件,就是指,我想提取的表型(本例子是小细胞肺癌SCLC)的所有样本的ID(什么是ID?求你看看GWAS分析-说人话(2)认识文件名)
然后,我们需要在表型数据(一般是一个txt文件上看,这个数据是怎样分的)
说人话:是的,右键打开来看就可以了!~
打开长这样:
然后,高逼格的就通过程序读##之后的第11行,低逼格的人工把前10行(有#的)去掉后保存成一个txt文件,如下:
本数据集中MORPHOLOGY==3,SCLC,所以我们就是提取所有第六行为3的ID出来。
好了,“大招”来了,打开你的R语言:
> tmp=read.table("phen.txt",header=T, stringsAsFactors = F, sep="\t")
> tmp[1:2,1:6]
> table(tmp$CASE_STATUS)
> table(tmp$MORPHOLOGY)
> smallcell=which(tmp$MORPHOLOGY==3)
> smallcell
> smallcellID=tmp$GENEVA_ID[smallcell]
> smallcellID
> write.table(smallcellID, "/Users/seedson/Desktop/SOM/phs000093.v2 data/SCLC_GENEVA_ID.txt", row.names= FALSE, col.names = FALSE)
这个操作就是把所有SCLC的ID全部提取出来,保存成一个SCLC_GENEVA_ID.txt的文件,没错,这个就下一步keep中间需要的文件了。
说人话:R代码解读如下:
> 读文件
> 看看有没有读取成功,第1-2行,1-6列长啥样
> 看看表格中CASE_STATUS
> 看看表格中MORPHOLOGY
> 定义smallcell是表格中MORPHOLOGY这一列等于3的数字
> 看看smallcell长啥样
> smallcellID就是表格中GENEVA_ID这一列中属于smallcell的数,提取出来
> 看看smallcellID长啥样
> 输出数据
我尽力了...(看得高兴就打个赏呗~)
上面的代码显示的整个分析过程,实际上每个组的代码就三行,如下:
下面就是重复的东西了:
接着,我还要对对照组,腺癌,鳞癌等的数据集提取出来(后面肯定是要分析不同的组是怎样变化的啦~):
按需求,反复多次:
鳞癌
> squamous_epi=which(tmp$MORPHOLOGY==1)
> squamous_epiID=tmp$GENEVA_ID[squamous_epi]
> write.table(squamous_epiID, "/Users/seedson/Desktop/SOM/phs000093.v2 data/squamous_epi_GENEVA_ID.txt", row.names= FALSE, col.names = FALSE)
腺癌
> adeno=which(tmp$MORPHOLOGY==2)
> adenoID=tmp$GENEVA_ID[adeno]
> write.table(adenoID, "/Users/seedson/Desktop/SOM/phs000093.v2 data/adeno_GENEVA_ID.txt", row.names= FALSE, col.names = FALSE)
大细胞
> large=which(tmp$MORPHOLOGY==4)
> largeID=tmp$GENEVA_ID[large]
> write.table(largeID, "/Users/seedson/Desktop/SOM/phs000093.v2 data/large_GENEVA_ID.txt", row.names= FALSE, col.names = FALSE)
数据集中其他的非小细胞肺癌
> otherNSCLC=which(tmp$MORPHOLOGY==7)
> otherNSCLCID=tmp$GENEVA_ID[otherNSCLC]
> write.table(otherNSCLCID, "/Users/seedson/Desktop/SOM/phs000093.v2 data/otherNSCLC_GENEVA_ID.txt", row.names= FALSE, col.names = FALSE)
数据中的其他肺癌
> otherlungca=which(tmp$MORPHOLOGY==8)
> otherlungcaID=tmp$GENEVA_ID[otherlungca]
> write.table(otherlungcaID, "/Users/seedson/Desktop/SOM/phs000093.v2 data/otherlungca_GENEVA_ID.txt", row.names= FALSE, col.names = FALSE)
对照组
> controlgroup=which(tmp$CASE_STATUS==0)
> controlgroupID=tmp$GENEVA_ID[controlgroup]
> write.table(controlgroupID, "/Users/seedson/Desktop/SOM/phs000093.v2 data/controlgroup_GENEVA_ID.txt", row.names= FALSE, col.names = FALSE)
最后,相应的文件夹(目录)中就会产生相应的数据: