======DeconSeq 简介======
从不纯核酸制剂中获得的序列可能包含来自样本以外来源的DNA。这些序列污染严重关注用于下游分析的数据质量，可能导致序列串联的错误组合和错误结论。因此，清除序列污染物是所有元基因组项目的必要步骤。

DeconSeq的交互式Web界面可用于自动检测和高效地从基因组和元基因组数据集中删除序列污染。[目前看来web暂时无法登陆]

可参考官方文档：
[[https://deconseq.sourceforge.net/manual.html#DB]]

[[https://sourceforge.net/projects/deconseq/files/]]

现整理了相关流程步骤如下：
======自制数据库构建======
下载序列数据
有几个地方可以检索序列数据。人类参考基因组构建37可以从国家生物技术信息中心（NCBI）FTP服务器下载。如果您收到错误，请检查文件在指定的位置是否仍然可用，以及文件名是否仍然相同。

  for i in {1..22} X Y MT; do wget ftp://ftp.ncbi.nih.gov/genomes/H_sapiens/Assembled_chromosomes/seq/hs_ref_GRCh37.p2_chr$i.fa.gz; done


提取和连接数据
序列数据使用GZIP算法和25个不同的文件进行压缩。下一个命令将提取数据并将其写入单个文件。同时，压缩文件将被删除（rm命令）。gzip命令中的-v提供了查看进度的简单选项。

  for i in {1..22} X Y MT; do gzip -dvc hs_ref_GRCh37.p2_chr$i.fa.gz >>hs_ref_GRCh37_p2.fa; rm hs_ref_GRCh37.p2_chr$i.fa.gz; done
----------------------------------------分割线----------------------------------------

对于物种的基因组而言，不需要执行上述的步骤，只需要把数据下载好，按照下面的方式进行处理即可

通过模糊基数N的长重复来分割序列
BWA和BWA-SW用随机基（ACGT）取代模糊基（所有不是ACGT的东西）。例如，这可能导致长时间的Ns中的假阳性命中，这些命中被随机转换为类似于查询序列的序列。替换模棱两可的基的原因是序列数据的2位表示，该表示只允许存储/表示四个不同的基（如00、01、10和11）。因此，在创建数据库之前，删除长长的模糊碱基是一个好主意（特别是考虑到人类参考基因组中50,000 Ns的延伸）。以下命令将额外删除每个序列开头和结尾的Ns（例如，1号染色体开头的10,000 Ns）。

  cat hs_ref_GRCh37_p2.fa | perl -p -e 's/N\n/N/' | perl -p -e 's/^N+//;s/N+$//;s/N{200,}/\n>split\n/' >hs_ref_GRCh37_p2_split.fa; rm hs_ref_GRCh37_p2.fa


过滤序列
BWA-SW在生成数据库时使用所有输入序列的串联。拆分基因组序列后，查询序列可能会比数据库中的序列更长。当查询序列在数据库中的多个（串联）序列上对齐时，这可能会导致问题。在当前版本中，BWA-SW无法解析跨越两个以上数据库序列的对齐。通过删除较短的序列，查询序列应该无法与数据库中的两个以上序列对齐。
除了短序列外，还应该删除重复的序列。重复不会添加新信息，但会增加数据库的大小。为了减少假阳性对齐的数量（见上文），应删除具有太多模糊基数的序列。
这些序列可以用PRINSEQ（http://prinseq.sourceforge.net）等程序轻松过滤。以下命令将额外重命名序列标识符，以确保整个数据集中的唯一标识符，并删除输入文件（rm命令）。

  perl /TJPROJ1/META_ASS/PreSaleEvaluation/remove_contaminations/deconseq-standalone-0.4.3/prinseq-lite.pl -log -verbose -fasta hs_ref_GRCh37_p2_split.fa -min_len 200 -ns_max_p 10 -derep 12345 -out_good hs_ref_GRCh37_p2_split_prinseq -seq_id hs_ref_GRCh37_p2_ -rm_header -out_bad null; rm hs_ref_GRCh37_p2_split.fa


创建数据库
BWA程序提供了一个选项，以FASTA格式从序列数据创建数据库。对于更大的数据集，应该使用bwtsw算法。数据库的最大大小为4GB。如果您想从更大的数据集生成数据库，则必须将数据拆分为小于4GB的块。索引过程可能需要一段时间，因此，我们将使用命令末尾的“&”在后台运行该过程。为了了解BWA在做什么，我们将将其输出写入bwa.log文件。您可以使用“$ tail -f bwa.log”来检查当前状态。

  /TJPROJ1/META_ASS/PreSaleEvaluation/remove_contaminations/deconseq-standalone-0.4.3/bwa64 index -p hs_ref_GRCh37_p2 -a bwtsw hs_ref_GRCh37_p2_split_prinseq.fasta >bwa.log 2>&1 &

当BWA程序完成后，您应该有八个由-p参数指定的名称的文件。

该步由于需要多次执行，建议批量生成，可以使用python3 deconseq_to_get_db.py > w.sh ，修改自己的路径或者增加参数后，配合 
  perl /PUBLIC/software/MICRO/share/MetaGenome_pipeline/MetaGenome_pipeline_V5.1/lib/00.Commbin/super_worker.pl --qalter --cyqt 1 --maxjob 200 --sleept 600  --qopts ' -q meta_ass.q ' w.sh --splits '\n\n' --prefix db --resource vf=1G superworker一起使用



======数据准备======
**分析前配置文件处理：**
需要修改/TJPROJ1/META_ASS/PreSaleEvaluation/remove_contaminations/deconseq-standalone-0.4.3/DeconSeqConfig.pm内的内容

{{:yuxi:企业微信20231010-222122_2x.png?400|}}

修改路径地址以及bwa所在的路径地址

增加数据库内的字典文件

{{:企业微信20231010-222306_2x.png?400|}} 

需要按照
  bact => {name => 'Bacterial genomes',
          db => 'bactDB'},
这样的格式添加数据库字典，其中第一个bact名称是后面deconseq.pl脚本需要用的到-dbs 参数，name不重要，db为使用bwa建立的索引名称
如果遇到多个数据库需要添加，可以使用路径下的get_db_dict.sh,直接生成字典，并复制到配置文件中
<code>
#!/bin/bash

names=("Azospirillum_aestuarii"
        "Brucella_sp."
        "Chryseobacterium_sp."
        "Cupriavidus_sp."
        "Trinickia_sp.")

for name in "${names[@]}"; do
    echo "'$name' => {name => '$name genomes',"
    echo "          db => '$name'},"
done
</code>


======流程======
配置好文件后，使用awk获得脚本命令
<code>
cd 目录
perl /TJPROJ1/META_ASS/PreSaleEvaluation/remove_contaminations/deconseq-standalone-0.4.3/deconseq.pl -f 需要剔除的序列 -dbs 物种数据库索引名称
</code>

获得文件为编号加clean或cont，其中clean是去除污染后的序列，cont是污染序列

====== 参考文献 ======
1. Li R, Li Y, Zheng H, Luo R, Zhu H, Li Q, Qian W, Ren Y, Tian G, Li J, Zhou G, Zhu X, Wu H, Qin J, Jin X, Li D, Cao H, Hu X, Blanche H, Cann H, Zhang X, Li S, Bolund L, Kristiansen K, Yang H, Wang J, Wang J: Building the sequence map of the human pan-genome. Nat. Biotechnol 2010, 28:57-63.

2. Li H, Durbin R: Fast and accurate long-read alignment with Burrows-Wheeler transform. Bioinformatics 2010, 26:589-595.