文章测序数据下载碰到的问题的小总结(南京站学员分享)

2019-10-25 00:46:19 浏览数 (1)

全国巡讲南京站过去还不到一周,学员们课后练习都很拼,微信群答疑对话一不留神就几百条了,爱学习的你是最美丽的!

恰好看到一个学员开始主动思考,自行摸索,超出我们授课范围的知识点整理,主动投稿,下面请大家欣赏一下南京站学员分享:

我目前的项目主要参考了一篇文献,我想借重复这篇文章的数据分析及表格制作来巩固学习3天的生信技能树培训内容。我翻到文章最后找到了数据的获取方式:

载入下载页面后一看有两个“数据大小”:Bytes和Bases,并且数据都不小,我电脑硬盘只有1T,如果下载数据是后者,那我连数据下载的第一步都无法完成,得先去扩展硬盘了。

那么这两个数字哪个是下载数据大小呢?

最直接的方式是将所有文件全部下载下来,看看文件的总大小,但是如果这样做耗时太长。我打算下载一个文件看看大小,然后找到tablelist里与文件大小对应列,统计一下整列的值的和,就得到了所有样本的数据量。看看这个数据量是不是对应的这两个数字中的一个。下载命令如下:

代码语言:javascript复制
nohup wget -c https://sra-download.ncbi.nlm.nih.gov/traces/sra76/SRR/008301/SRR8500494 &
ls -lh

我下载的样本名称编号是SRR8500494,这个样本文件的大小是549M

在页面下方的表格中找到对应的样本,该数字对应的列名是MBytes,在R中统计MBytes和MBases这两列的数据的总和,

命令及结果如下:

代码语言:javascript复制
> rm(list = ls())
> options(stringsAsFactors = F)
> x <- read.table("SraRunTable.txt",sep = ',',header = T)
> dim(x)
[1] 119  33
> colnames(x)
 [1] "Run"                  "Assay.Type"          
 [3] "AvgSpotLen"           "Biomaterial_provider"
 [5] "BioProject"           "BioSample"           
 [7] "BioSampleModel"       "Center.Name"         
 [9] "Consent"              "DATASTORE.filetype"  
[11] "DATASTORE.provider"   "DATASTORE.region"    
[13] "Experiment"           "Instrument"          
[15] "isolate"              "Library.Name"        
[17] "LibraryLayout"        "LibrarySelection"    
[19] "LibrarySource"        "MBases"              
[21] "MBytes"               "Organism"            
[23] "Platform"             "ReleaseDate"         
[25] "sample_acc"           "Sample.Name"         
[27] "sex"                  "SRA.Study"           
[29] "tissue"               "disease"             
[31] "Age"                  "capture"             
[33] "replicate"           
> y <- x$MBytes
> z <- x$MBases
> sum(y)
[1] 314988
> sum(z)
[1] 915887

结果说明Bytes的数值代表下载下来的SRA文件的大小。

(最开始并不知道页面里可以勾选单个样本,并且页面会显示样本的这2个数据)

插曲1:刚开始没有加参数,结果报错如下

代码语言:javascript复制
> rm(list = ls())
> options(stringsAsFactors = F)
> x <- read.table("SraRunTable.txt")
Error in scan(file = file, what = what, sep = sep, quote = quote, dec = dec,  : 
  line 1 did not have 14 elements

后来用Notepad 打开,发现数据分隔符是“,”,猜想可能是分隔符的问题,搜索后添加打开函数的参数,结果成功打开。

插曲2:用R语言统计的Bytes的结果是314988M,与307.61 Gb看上去不相等,实际是以1024为倍数换算是一致的;用R语言统计的Bases的结果是915887,与915.89 G确是以1000换算的。

再看看shell里的换算倍数是哪个:

嗯,是1024。

另外一个值的单位是bp,我猜想是数据读到的碱基数总和,为了验证我的猜想,我将下载下来的文件进行了SRA->Fastq文件的转换,并想用Fastqc统计总base数,结果发现Fastqc并不会统计这个值,于是在网上搜索了如下一行命令直接可以得到这个值:

代码语言:javascript复制
(base) $ cat SRR8500494_1.fastq | paste - - - - | cut -f 4 | fold -w 1 | wc -l
882574875
(base) $ cat SRR8500494_2.fastq | paste - - - - | cut -f 4 | fold -w 1 | wc -l
882574875

可以看到,R1和R2的值相加刚好是该样本对应的Bases值:1.76G。

至此,这个两个数值是代表什么就搞清楚了:Bytes代表下载下来的SRA文件的大小,Bases代表的是这个数据读取的碱基数总和

不知道你有没有发现,就是没有人提到fastq文件的大小呢?

聪明的读者朋友,可以留言说出自己的看法,为什么不标记fastq文件的大小?

另外,推荐学徒数据挖掘第二期汇总之多分组基因注释代码大放送(compareCluster) ,部分精彩推文如下:

  • 热图、韦恩图、GO富集分析图
  • 纯R代码实现ssGSEA算法评估肿瘤免疫浸润程度
  • 按基因在染色体上的顺序画差异甲基化热图
  • 有生物学意义的复杂热图
  • 干扰MYC-WWP1通路重新激活PTEN的抑癌活性——3步搞定GSEA分析
  • TCGA正常血液样本中PD-L1与BRCA1和NBS1的表达量相关性
  • 肿瘤异质性 免疫浸润细胞数据挖掘(可能是最简单的3分文章了)

0 人点赞