比对其实就是如此
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
sam/bam/pos
- 上述identifier由比对结果看,位于chr1;
- 从fasta中将chr1序列提取出来,截取出其对应的碱基位置,进行查看;
- 依据,SAM中对于pos的解释,为比对上的reference左端的位置,所以,最后截取的reference的位置,为
左(小)~(左(大)+读长-1)
,这也就是echo $((110954876+149))
的来由;grep -i
所跟序列为SAM|BAM中所对应的比对序列;
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