update-20250601

0_raw/step.sh

@@ -0,0 +1,9 @@
+###从NCBI下载基因组(DNA)/转录组(RNA)测序数据到获得预测基因的处理步骤（以Nyctotherus_ovalis物种为例，见/home/Data_01/yizhenzhen06/public_share/Saccharomyces_cerevisiae_DNA。如果没有注明仅基因组数据要做，其他步骤基因组和转录组步骤一样。直接写在本文件里的美丽表示可以直接终端运行指令，其他指令要提交pbs任务运行。
+
+screen -S ncbi_dl #screen命令创建一个新的窗口（这里窗口命名为ncbi_dl）用来下载原始测序数据
+wget -c -b -i species_name_urls.txt #NCBI下载每个物种相应所有原始数据至服务器，species_name_urls.txt包括SRA数据的urls链接
+#命令运行成功后，按下组合键ctrl+A+D退出当前窗口
+screen -r ncbi_dl #返回ncbi_dl窗口，查看下载进度
+
+:fastq-dump --defline-seq '@$sn[_$rn]/$ri' --split-files SRR10355986.1 #RNA数据转换数据格式为左右两端read文件
+:fastq-dump --split-files SRR10355986.1 #DNA数据转换数据格式为左右两端read文件

1_trim/step.sh

@@ -0,0 +1,23 @@
+###从NCBI下载基因组(DNA)/转录组(RNA)测序数据到获得预测基因的处理步骤（以Nyctotherus_ovalis物种为例，见/home/Data_01/yizhenzhen06/public_share/Saccharomyces_cerevisiae_DNA。如果没有注明仅基因组数据要做，其他步骤基因组和转录组步骤一样。直接写在本文件里的美丽表示可以直接终端运行指令，其他指令要提交pbs任务运行。
+
+fastqc -o 1_trim/QC/ --noextract -t 2 -f fastq 0_rawdata/SRRXXXXXX_1.fastq 0_rawdata/SRRXXXXXX_2.fastq
+
+#分别检查左、右端数据质量，数据量大可以一起提交脚本运行。
+#如果数据质量不好，就做trimmomatic质控过滤。如果质量没问题，就分别合并所有左、右端两端reads文件进行组装。
+
+:cat SRRXXXXXX.1_1.fastq SRRXXXXXX.1_1.fastq > species_name_1.fastq #合并所有左端reads文件
+:cat SRRXXXXXX.1_2.fastq SRRXXXXXX.1_2.fastq > species_name_2.fastq #合并所有右端reads文件
+
+trimmomatic PE -phred33 -threads 16 \
+1_trim/species_name_1.fastq \
+1_trim/species_name_2.fastq \
+1_trim/species_name_1_p.fastq \
+1_trim/species_name_1_up.fastq \
+1_trim/species_name_2_p.fastq \
+1_trim/species_name_2_up.fastq \
+ILLUMINACLIP:../../adapters/NexteraPE-PE.fa:2:30:10:8:true LEADING:3 TRAILING:3 SLIDINGWINDOW:5:20 MINLEN:36 HEADCROP:15 
+#trimmomatic质控过滤，参数需根据具体数据质量进行调整，大型数据请提交脚本运行。
+
+fastqc -o 1_trim/QC/ --noextract -t 2 -f fastq 0_rawdata/SRRXXXXXX_1.fastq 0_rawdata/SRRXXXXXX_2.fastq
+
+#再次检查左、右端数据质量，数据量大同样提交脚本运行。

2_assembly/step.sh

@@ -0,0 +1,18 @@
+###从NCBI下载基因组(DNA)/转录组(RNA)测序数据到获得预测基因的处理步骤（以Nyctotherus_ovalis物种为例，见/home/Data_01/yizhenzhen06/public_share/Saccharomyces_cerevisiae_DNA。如果没有注明仅基因组数据要做，其他步骤基因组和转录组步骤一样。直接写在本文件里的美丽表示可以直接终端运行指令，其他指令要提交pbs任务运行。
+
+Trinity --seqType fq --max_memory 50G --no_version_check \
+--left species_name_1_p.fastq \
+--right species_name_2_p.fastq \
+--CPU XX 
+#RNA数据用Trinity组装，--CPU调整调用cpu核心数，其它参数需根据具体数据情况进行调整，提交脚本运行。
+
+spades.py --sc -t XX -m 230 \
+-1 species_name_1_p.fastq \
+-2 species_name_2_p.fastq \
+-o spades_species_name
+#DNA数据用spades组装，-t调整调用cpu核心数，其它参数需根据具体数据情况进行调整，提交脚本运行。
+
+busco -i spades/scaffolds.fasta -m XXX \
+-l /home/Data_01/share_databases/busco/alveolata_odb10/alveolata_odb10 \
+-o busco_out_alve --offline -f 
+#BUSCO 检测完整度，-m调整评估模式(trans转录组 geno基因组)

3_decontamination/step.sh

@@ -0,0 +1,58 @@
+###从NCBI下载基因组(DNA)/转录组(RNA)测序数据到获得预测基因的处理步骤（以Nyctotherus_ovalis物种为例，见/home/Data_01/yizhenzhen06/public_share/Saccharomyces_cerevisiae_DNA。如果没有注明仅基因组数据要做，其他步骤基因组和转录组步骤一样。直接写在本文件里的美丽表示可以直接终端运行指令，其他指令要提交pbs任务运行。
+
+blastn -db /home/Data_01/share_databases/ncbi_genome_refseq_bac_old/ncbi_genome_refseq_bac_old \
+-query /home/Data_01/yizhenzhen06/public_share/Saccharomyces_cerevisiae_DNA/2_assembly/spades/scaffolds.fasta \
+-out 1_Sacc_cerevisiae_dec_bac.tab \
+-outfmt "6 qseqid qlen sseqid slen pident qcovs qstart qend sstart send mismatch gapopen evalue bitscore" \
+-evalue 1e-5 -max_target_seqs 1 -num_threads X
+
+sed -i '1i\tquery ID\tquery sequence length\tsubject ID\tsubject qequence length\tpercent identities\tcoverage\tqstart\tqend\tsstart\tsend\tmismatch\tgapopen\tevalue\tbitscores ' 1_Sacc_cerevisiae_dec_bac.tab
+#去细菌基因组污染，两条指令放一个脚本运行（注意检查数据库路径）
+
+awk '{ if ($5 > 80 && $6 > 70) print $1 }' 1_Sacc_cerevisiae_dec_bac.tab > 1_Sacc_cerevisiae_dec_bac.id #提取第5列相似度大于80%且覆盖率大于70%的序列ID（这一步建议需要下载到本地手动检查输出的id条数是否和Sacc_cerevisiae_dec_bac.tab（手动excel筛选）文件的第5列相似度大于80%且覆盖率大于70%的序列条数一致）
+seqkit grep -v -f Sacc_cerevisiae_dec_bac.id ../2_assembly/trinity_out_dir/Trinity.fasta > 1_Sacc_cerevisiae_dec_bac.fasta #从组装结果文件中删去细菌污染序列
+
+####################################################################################################
+
+blastn -db /home/Data_01/share_databases/human_genome_new/GCA_021234545.1_HCC1395BL_v1.0_genomic_db \
+-query 1_Sacc_cerevisiae_dec_bac.fasta \
+-out 2_Sacc_cerevisiae_dec_bac_hum.tab \
+-outfmt "6 qseqid qlen sseqid slen pident qcovs qstart qend sstart send mismatch gapopen evalue bitscore" \
+-evalue 1e-5 -max_target_seqs 1 -num_threads X
+
+sed -i '1i\tquery ID\tquery sequence length\tsubject ID\tsubject qequence length\tpercent identities\tcoverage\tqstart\tqend\tsstart\tsend\tmismatch\tgapopen\tevalue\tbitscores ' 2_Sacc_cerevisiae_dec_bac_hum.tab
+#去人类基因组污染，两条指令放一个脚本运行（注意检查数据库路径）
+
+awk '{ if ($5 > 80 && $6 > 70) print $1 }' 2_Sacc_cerevisiae_dec_bac_hum.tab > 2_Sacc_cerevisiae_dec_bac_hum.id #提取第5列相似度大于80%且覆盖率大于70%的序列ID（这一步建议需要下载到本地手动检查输出的id条数是否和2_Sacc_cerevisiae_dec_bac_hum.tab（手动excel筛选）文件的第5列相似度大于80%且覆盖率大于70%的序列条数一致）
+seqkit grep -v -f 2_Sacc_cerevisiae_dec_bac_hum.id 1_Sacc_cerevisiae_dec_bac.fasta > 2_Sacc_cerevisiae_dec_bac_hum.fasta #从组装结果文件中删去人类基因组污染序列
+
+####################################################################################################
+
+blastn -db /home/Data_01/share_databases/human_genome_new/GCF_000001405.39_GRCh38.p13_rna_db \
+-query 2_Sacc_cerevisiae_dec_bac_hum.fasta \
+-out 3_Sacc_cerevisiae_dec_bac_hum_geno_rna.tab \
+-outfmt "6 qseqid qlen sseqid slen pident qcovs qstart qend sstart send mismatch gapopen evalue bitscore" \
+-evalue 1e-5 -max_target_seqs 1 -num_threads X
+
+sed -i '1i\tquery ID\tquery sequence length\tsubject ID\tsubject qequence length\tpercent identities\tcoverage\tqstart\tqend\tsstart\tsend\tmismatch\tgapopen\tevalue\tbitscores ' 3_Sacc_cerevisiae_dec_bac_hum_geno_rna.tab
+#去人类转录组污染，两条指令放一个脚本运行（注意检查数据库路径）
+
+awk '{ if ($5 > 80 && $6 > 70) print $1 }' 3_Sacc_cerevisiae_dec_bac_hum_geno_rna.tab > 3_Sacc_cerevisiae_dec_bac_hum_geno_rna.id #提取第5列相似度大于80%且覆盖率大于70%的序列ID（这一步需要下载到本地手动检查输出的id条数是否和Nyct_bac_human_genome_blastn.tab（excel筛选）文件的第5列相似度大于80%且覆盖率大于70%的序列条数一致）
+seqkit grep -v -f 3_Sacc_cerevisiae_dec_bac_hum_geno_rna.id 2_Sacc_cerevisiae_dec_bac_hum.fasta > 3_Sacc_cerevisiae_dec_bac_hum_geno_rna.fasta #从组装结果文件中删去人类转录组污染序列
+
+####################################################################################################
+
+seqkit seq -n -i ../2_assembly/spades/scaffolds.fasta > 4_Sacc_cerevisiae_dec_all.txt
+#提取序列ID，仅基因组数据要做
+awk '{if($1~"_0.")print$1}' 4_Sacc_cerevisiae_dec_all.txt > 4_Sacc_cerevisiae_dec_all_cov1.id
+#提取coverage小于1的序列ID，仅基因组数据要做
+seqkit grep -v -f 4_Sacc_cerevisiae_dec_all_cov1.id 4_Sacc_cerevisiae_dec_all.fasta > 4_Sacc_cerevisiae_dec_all_cov1.fasta 
+#去除coverage小于1的序列，仅基因组数据要做
+
+####################################################################################################
+
+cd-hit-est \
+-i 4_Sacc_cerevisiae_dec_all_cov1.fasta \
+-o 5_Sacc_cerevisiae_dec_all_cov1_cdhit.fasta \
+-c 0.9 -n 8 -d 0 -M 32000 -T 48
+#去冗余（核酸序列）

4_prediction/step.sh

@@ -0,0 +1,34 @@
+###从NCBI下载基因组(DNA)/转录组(RNA)测序数据到获得预测基因的处理步骤（以Nyctotherus_ovalis物种为例，见/home/Data_01/yizhenzhen06/public_share/Saccharomyces_cerevisiae_DNA。如果没有注明仅基因组数据要做，其他步骤基因组和转录组步骤一样。直接写在本文件里的美丽表示可以直接终端运行指令，其他指令要提交pbs任务运行。
+
+seqkit sample ../3_decontamination/5_Sacc_cerevisiae_dec_all_cov1_cdhit.fasta \
+-o 5_Sacc_cerevisiae_dec_all_cov1_cdhit_sample.fasta \
+-p 0.25 -s 11 
+#seqkit抽提 1/2 或 1/4 序列 -p 参数调整
+
+#Codetta 读取基因组，然后利用已知蛋白质的数据库，通过比对保守蛋白质序列的隐马尔可夫模型（Profile HMMs）与输入序列，计算出每个密码子最可能对应的氨基酸 
+#如果一个基因编码某种蛋白质，该蛋白质的氨基酸序列应该与已知蛋白质具有一定的相似性，Codetta 就通过寻找这种相似性来推断每个密码子的氨基酸解码
+
+codetta.py 5_Sacc_cerevisiae_dec_all_cov1_cdhit_sample.fasta 
+#codetta预测密码子表（进入 py36 环境，提交任务）
+
+####################################################################################################
+
+augustus \
+--progress=true --protein=on --introns=on --start=on --stop=on --cds=on --codingseq=on --gff3=on \
+--species=tetrahymena \
+--extrinsicCfgFile=/home/Data_01/yizhenzhen06/miniconda3/envs/augustus/config/extrinsic/extrinsic.M.RM.E.W.cfg \
+../3_decontamination/5_Sacc_cerevisiae_dec_all_cov1_cdhit.fasta \
+--outfile=Sacc_cerevisiae.gff3
+#基因组数据基因预测，提交脚本运行
+getAnnoFasta.pl Sacc_cerevisiae.gff3 && mv Sacc_cerevisiae.aa Sacc_cerevisiae_augustus.fasta
+# 从预测结果中提取氨基酸序列
+
+TransDecoder.LongOrfs \
+-t ../3_decontamination/5_Sacc_cerevisiae_dec_all_cov1_cdhit.fasta \
+-O LongOrfs/Sacc_cerevisiae_clean_cdhit_ciliate \
+-G Ciliate
+TransDecoder.Predict \
+-t ../3_decontamination/5_Sacc_cerevisiae_dec_all_cov1_cdhit.fasta \
+-O Sacc_cerevisiae_clean_cdhit_predict \
+-G Ciliate
+#转录组数据基因预测，-G参数改成help命令中相应物种的对应密码子，两条指令一起提交脚本运行