bilibili视频-《转录组测序数据分析》笔记

接下来要从fa转fq文件

ls -lh raw/|wc #看数据
less raw/nohup.out #有后台日志
grep failed raw/nohup.out
zless SRR1039508_1.fastq.gz #查看

非常长,要用报告来看,如何做报告,如下:

fastqc SRR1039510_1.fastq.gz
fastqc -t -10 SRR1039510_1.fastq.gz#-t -10加参数使速度加快

(写不进去,是因为是老大的命令,要放进自己的路径)

ls *gz |xargs fastqc -t 10#批量生成fastqc报告

fastqc文件生产后,用multiqc,综合生产一个报告

multiqc ./

老大视频里的未修改前的循环

更改后的

在后台运行

注:如果有config,在cat后改为config,如果没有就是$1

source activate rna
bin_trim_galore=trim_galore
dir='/mm/project/airway/clean'
cat $1 |while read id
do
      arr=($id)
      fq1=${arr[0]}
      fq2=${arr[1]}
nohup $bin_trim_galore -q 25 --phred33 --length 36 -e 0.1 --stringency 3 --paired -o $dir $fq1 $fq2 &
done
source deactivate
dir='/mm/project/airway/clean'
cat config |while read id
do
      arr=($id)
      fq1=${arr[0]}
      fq2=${arr[1]}
nohup trim_galore -q 25 --phred33 --length 36 -e 0.1 --stringency 3 --paired -o $dir $fq1 $fq2 &
done

trim_galore一个

之前就没跟着一起跑过这个循环,自己弄也跑不出来,和老大的代码一样也不行,改下也不行,真不知道怎么回事了.....就是没有我的任务😢

config也没啥问题

划重点了,我要默念一百遍!以上是我未跑出来的过程,为了记录过程,还是保留。最后还是我 老大一语中的,我明明是就改了dir的输出路径,出了问题,那我就应该把问题集中在路径这里就好了!

###这个是最后能跑的

source activate rna
bin_trim_galore=trim_galore
dir='/four/mm/project/airway/clean'  ###问题就在这个路径这,我没从根目录出发,是从mm这个相对路径出发的
cat $1 |while read id
do
      arr=($id)
      fq1=${arr[0]}
      fq2=${arr[1]}
nohup $bin_trim_galore -q 25 --phred33 --length 36 -e 0.1 --stringency 3 --paired -o $dir $fq1 $fq2 &
done
source deactivate

已经跑了

好像又出来什么问题,任务显示能跑,但是为什么是gzip

单个跑

trim_galore --phred33 -q 25 -e 0.1 --length 36 --stringency 3 --paired -o ./ ~/fqmm/SRR1039508_1.fastq.gz ~/fqmm/SRR1039508_2.fastq.gz

RNA-seq:7-alignment

从老大视频里扒来的链接,先收藏

https://blog.csdn.net/xubo245/article/details/50878760

https://blog.csdn.net/xubo245/article/details/50836185

比对:hisat2、subjunc、star大多是针对转录组的,bwa、bowtie2是基因组

老大的参考基因组位置

/public/reference/genome/*

老大的参考基因组索引位置

/Public/reference/index/*

加了✳️和不加是不一样的,加了✳️可以把文件夹内的内容同时显示出来,不加的话只会显示当前路径下的文件夹

hisat有单独文件夹

/public/reference/index/hisat/*

每个文件取前1000行

ls ./*gz |while read ijfmklmf;do(zcat $ijfmklmf |head -1000 > 一个文件名);done

因为加了管道符,所以要加()

由于输出到clean文件里了,所以想只要文件名

改名字

sam文件

我跟着做的

是不是链接过来的文件不能做下面这样的操作呢?

上面是链接过来的文件就不行,我自己trim_galore的文件就可以,不知道为什么啊

上面是链接过来的文件就不行,我自己trim_galore的文件就可以

但是改成去掉“gz”的以后,就不能比对了😢哇,好神奇,知道了呢!是因为我按照老大的去掉了“.gz”后,就把fq.gz的文件给移动到别的文件夹了,然后就不可以了,因为RR1039508_1_val_1.fq只是SRR1039508_1_val_1.fq.gz的类似快捷方式的东西吧,并不是文件本身,-lh也能看到它没有大小。我在把fq.gz文件移动回来以后,就可以啦

呃,但是报错了,是因为fq-2是0,这是为什么呢/😢

视频里老大的比对:

单个文件比对:

hisat2 -p 4 -x /four/mm/index/hisat/hg38/genome -1 SRR1039508_1_val_1.fq.gz -2 SRR1039508_2_val_2.fq.gz -S SRR1039508.tmp.sam#也可以输出为SRR1039508.hisat.sam,是没去掉gz的,要去掉gz的原因是有些软件会识别成压缩软件 

总之,不去掉.gz是可以比对的,就行了😄

用循环比对,hisat2比对到索引文件,然后samtools直接对生成的bam文件排序,以利于后续软件分析。

nohup cat SRR_Acc_List.txt |while read id;do  #复制一份id到当前路径下
hisat2 -p 5 -x ~/index/grch38/genome -1 \
${id}_1_val_1.fq.gz -2 ${id}_2_val_2.fq.gz | \
samtools sort -@ 5 -o ~/rna.GSE52778/sort.bam/${id}.sort.bam -
done &

差异分析

崔老师ppt里的,也没跑出来

for fn in {508..523}
do
featureCounts -T 5 -p -t exon -g gene_id \
-a /four/mm/project/gtf/gencode.v29.annotation.gtf \ 
-o SRR1039$fn.counts.txt SRR1039$fn.sort.bam
# donot set dir
done

###我改完以后的,跑不出来
cat SRR_Acc_List.txt |while read id;do featureCounts -T 5 -p -t exon -g gene_id 
-a /four/mm/project/gtf/gencode.v29.annotation.gtf 
-o ${id}.counts.txt ${id}.sort.bam 
# donot set dir 
done
最后编辑于
©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 216,470评论 6 501
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 92,393评论 3 392
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 162,577评论 0 353
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 58,176评论 1 292
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 67,189评论 6 388
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 51,155评论 1 299
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 40,041评论 3 418
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 38,903评论 0 274
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 45,319评论 1 310
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 37,539评论 2 332
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 39,703评论 1 348
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 35,417评论 5 343
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 41,013评论 3 325
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 31,664评论 0 22
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 32,818评论 1 269
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 47,711评论 2 368
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 44,601评论 2 353

推荐阅读更多精彩内容

  • 生信学习笔记 linux部分功能 查看文件夹 工具 选项 可以设置鼠标功能 可以设置右键粘贴 双击这个窗口可以再打...
    Vikenn阅读 1,150评论 1 4
  • 继续装一下虚拟环境: 因为有beautifulsoup4所以出现了deprecation是正常的。书中也提到了。 ...
    每日派森阅读 234评论 0 0
  • 项目开发中有的时候第三库需要版本号高的Ruby的版本,需要指定不同的版本. 查看所有的Ruby版本 2.查看当前R...
    FlyElephant阅读 2,694评论 0 3
  • 2018年3月23日 星期五 晴 今天下班回家闺女带着大红花回来了,一问才知道在学...
    张偌鸿亲子日记阅读 86评论 0 0
  • 心安,则一切皆安。 这不是一个让人心安的时代。 生存竞争激烈,房价高企,阶层固化,故土丧失,环境污染,实利至上。 ...
    知识管理某李阅读 1,000评论 8 2