老师,您好,根据咱们转录组课程进行了实操,质控等正常,但是转录组数据比对到基因组这个步骤总是出现比对率只有零点几的情况,尝试了好几个转录组数据,从来没正常过。
水稻转录组测序数据显示如下,NCBI下载的单端数据:1 @SRR22419667.1 A00133:246:HNTYVDSXY:1:1101:30734:1454 length=300
2 TCCCGTTGAGGTTACACTCCAAATTGAATGATGAGGCAACGAATAACAAAGCAGCACTCAGAAATCCATGATGGTCCAGGCATCTACAAATGATATATCGCGCAACTCCATTTGAGTACAAATGGGTCCAAGCCCTGGCAACAGGTGTATCCTTAAGAG
3 +SRR22419667.1 A00133:246:HNTYVDSXY:1:1101:30734:1454 length=300
4 FFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFF
5 @SRR22419667.2 A00133:246:HNTYVDSXY:1:1101:6822:1689 length=300
6 GAGCTCAGTAACAGCAGTCTCTGAGTCAGCATGCCATGTGAGTAAATCAACAAAAGGACAAGTACCAGGATCTAGCAAATTGGGTCAACTAGCTACACGGAGAAATTTAGAGACTGAAGTTTTGCAAAGTCAGCAAGTTAGCAGCCCAGCCAGTTCCAC
7 +SRR22419667.2 A00133:246:HNTYVDSXY:1:1101:6822:1689 length=300
8 FFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFF.
#水稻基因组的文件从ensemble plants上下载,构建索引的代码为
#基因组信息下载到服务器/home/us024/Imranrnaseq/my_rnaseq/ref
wget https://ftp.ensemblgenomes.ebi.ac.uk/pub/plants/release-61/fasta/oryza_sativa/dna/Oryza_sativa.IRGSP-1.0.dna.toplevel.fa.gz
wget https://ftp.ensemblgenomes.ebi.ac.uk/pub/plants/release-61/gff3/oryza_sativa/Oryza_sativa.IRGSP-1.0.61.gff3.gz
gunzip -c Oryza_sativa.IRGSP-1.0.dna.toplevel.fa.gz > Oryza_sativa.IRGSP-1.0.dna.toplevel.fa
gunzip -c Oryza_sativa.IRGSP-1.0.61.gff3.gz > Oryza_sativa.IRGSP-1.0.61.gff3
#修改gff3文件
sed 's#gene:##' Oryza_sativa.IRGSP-1.0.61.gff3|ee
sed 's#gene:##' Oryza_sativa.IRGSP-1.0.61.gff3|sed 's#transcript:##' |ee
#然后把修改后的另存为另一个名字,运行如下
sed 's#gene:##' Oryza_sativa.IRGSP-1.0.61.gff3|sed 's#transcript:##' > Oryza_sativa.IRGSP-1.0.61.gff31
#把源文件覆盖掉
mv Oryza_sativa.IRGSP-1.0.61.gff31 Oryza_sativa.IRGSP-1.0.61.gff3
#构建索引
sh $scriptdir/index.sh Oryza_sativa.IRGSP-1.0.dna.toplevel.fa Oryza_sativa.IRGSP-1.0.61.gff3
#构建索引后,gtf的部分显示,其它的生成的结果看着也没什么问题
1 1 RAP2022-09-01 transcript 2983 10815 . + . transcript_id "Os01t0100100-01"; gene_id "Os01g0100100";
2 1 RAP2022-09-01 exon 2983 3268 . + . transcript_id "Os01t0100100-01"; gene_id "Os01g0100100";
3 1 RAP2022-09-01 exon 3354 3616 . + . transcript_id "Os01t0100100-01"; gene_id "Os01g0100100";
4 1 RAP2022-09-01 exon 4357 4455 . + . transcript_id "Os01t0100100-01"; gene_id "Os01g0100100";
5 1 RAP2022-09-01 exon 5457 5560 . + . transcript_id "Os01t0100100-01"; gene_id "Os01g0100100";
6 1 RAP2022-09-01 exon 7136 7944 . + . transcript_id "Os01t0100100-01"; gene_id "Os01g0100100";
7 1 RAP2022-09-01 exon 8028 8150 . + . transcript_id "Os01t0100100-01"; gene_id "Os01g0100100";
8 1 RAP2022-09-01 exon 8232 8320 . + . transcript_id "Os01t0100100-01"; gene_id "Os01g0100100";
9 1 RAP2022-09-01 exon 8408 8608 . + . transcript_id "Os01t0100100-01"; gene_id "Os01g0100100";
10 1 RAP2022-09-01 exon 9210 9615 . + . transcript_id "Os01t0100100-01"; gene_id "Os01g0100100";
11 1 RAP2022-09-01 exon 10102 10187 . + . transcript_id "Os01t0100100-01"; gene_id "Os01g0100100";
12 1 RAP2022-09-01 exon 10274 10430 . + . transcript_id "Os01t0100100-01"; gene_id "Os01g0100100";
13 1 RAP2022-09-01 exon 10504 10815 . + . transcript_id "Os01t0100100-01"; gene_id "Os01g0100100";
14 1 RAP2022-09-01 CDS 3449 3616 . + 0 transcript_id "Os01t0100100-01"; gene_id "Os01g0100100";
15 1 RAP2022-09-01 CDS 4357 4455 . + 0 transcript_id "Os01t0100100-01"; gene_id "Os01g0100100";
16 1 RAP2022-09-01 CDS 5457 5560 . + 0 transcript_id "Os01t0100100-01"; gene_id "Os01g0100100";
#构建索引后,gff3的部分显示
69 #!genebuild-last-updated 2019-06
70 1 IRGSP-1.0 chromosome 1 43270923 . . . ID=chromosome:1;Alias=Chr1,AP014957.1,NC_029256.1
71 ###
72 1 RAP2022-09-01 gene 2983 10815 . + . ID=Os01g0100100;biotype=protein_coding;gene_id=Os01g0100100;logic_name=rapdb_genes
73 1 RAP2022-09-01 mRNA 2983 10815 . + . ID=Os01t0100100-01;Parent=Os01g0100100;biotype=protein_coding;tag=Ensembl_canonical;tra
74 1 RAP2022-09-01 exon 2983 3268 . + . Parent=Os01t0100100-01;Name=Os01t0100100-01-E1;constitutive=1;ensembl_end_phase=-1;ense
75 1 RAP2022-09-01 five_prime_UTR 2983 3268 . + . Parent=Os01t0100100-01
76 1 RAP2022-09-01 five_prime_UTR 3354 3448 . + . Parent=Os01t0100100-01
77 1 RAP2022-09-01 exon 3354 3616 . + . Parent=Os01t0100100-01;Name=Os01t0100100-01-E2;constitutive=1;ensembl_end_phase=0;ensem
78 1 RAP2022-09-01 CDS 3449 3616 . + 0 ID=CDS:Os01t0100100-01;Parent=Os01t0100100-01;protein_id=Os01t0100100-01
79 1 RAP2022-09-01 exon 4357 4455 . + . Parent=Os01t0100100-01;Name=Os01t0100100-01-E3;constitutive=1;ensembl_end_phase=0;ensem
80 1 RAP2022-09-01 CDS 4357 4455 . + 0 ID=CDS:Os01t0100100-01;Parent=Os01t0100100-01;protein_id=Os01t0100100-01
81 1 RAP2022-09-01 exon 5457 5560 . + . Parent=Os01t0100100-01;Name=Os01t0100100-01-E4;constitutive=1;ensembl_end_phase=2;ensem
#比对到基因组,先比对一个
hisat2 -p 3 --rg-id=400ppm_rep1 --rg SM:400ppm_rep1 --rg LB:400ppm_rep1 --rg PL:ILLUMINA -x $REF_INDEX --dta -U $workdir/2.data_qc/400ppm_rep1.clean.fq.gz -S 400ppm_rep1.sam 2> 400ppm_rep1.summary
#比对结果
24819002 reads; of these:
24819002 (100.00%) were unpaired; of these:
24775361 (99.82%) aligned 0 times
39911 (0.16%) aligned exactly 1 time
3730 (0.02%) aligned >1 times
0.18% overall alignment rate