生物信息学转录组shell命令linux

比对其实就是如此

2019-10-29  本文已影响0人  Juan_NF

fastq1和fastq2特定的cluster的序列

  • @SRR8980088.3 fastq的identifier
  • 基于此可以找到下机数据中该编号对应的序列;
a=/trainee/vip13/test
####查看对应的序列
zcat /trainee/vip13/test/SRR8980088_2_val_2.fq.gz|grep -A1 'SRR8980088.3'|head -n2|grep -v @
zcat /trainee/vip13/test/SRR8980088_1_val_1.fq.gz|grep -A1 'SRR8980088.3'|head -n2|grep -v @
fastq1 fastq2

fastq在SAM|BAM文件中的展示

  • 基于此,可以在SAM/BAM文件中找到,其比对到fasta上的序列
  • 在这里的体现是:fastq1与fa方向一致,fastq2与fa反向互补;
cat /trainee/vip13/test/SRR8980088.sam|grep SRR8980088.3|head -n2|cut -f1,10
SAM/BAM

Fastq中的read在fasta中的位置-SAM/BAM

  • 上述identifier由比对结果看,位于chr1;
  • 从fasta中将chr1序列提取出来,截取出其对应的碱基位置,进行查看;
  • 依据,SAM中对于pos的解释,为比对上的reference左端的位置,所以,最后截取的reference的位置,为 左(小)~(左(大)+读长-1),这也就是echo $((110954876+149))的来由;
  • grep -i所跟序列为SAM|BAM中所对应的比对序列;
sam/bam/pos
lib_PEseq.jpg
##从fa文件将chr1序列提取出来
zcat /trainee/vip13/test/hg38.fa.gz |sed -n '/chr1$/,/chr10$/p' > test.txt
###需要精确确定位置,所以,chr的字样需要去除
cat test.txt |grep -v '>chr' > chr1.txt
echo $((110954876+149))
###相应的,换行符等也要去掉
cat /trainee/vip13/test/chr1.txt |tr -d '\n'|cut -c 110954736-110955025|grep -i 
sam_fq1
sam_fq2
上一篇下一篇

猜你喜欢

热点阅读