寫在前面

在進行染色體掛載之前一般對得到的 "primary assembly"（主要組裝/初始組裝結(jié)果）進行進一步的優(yōu)化，以減少錯誤，提高基因組的質(zhì)量和準(zhǔn)確性绕沈。然后再將這些相對碎片化但polish過的高質(zhì)量primary assembly進行染色體掛載以得到一個染色體水平的組裝結(jié)果道偷。如果用CLR數(shù)據(jù)組裝出的primary assembly在沒有進行polish的情況下廷臼，直接用HIC數(shù)據(jù)進行染色體掛載會有如下問題：（以下幾點回答來自gpt3.5）

1. 錯誤和噪聲： CLR（Continuous Long Reads寡痰，如 PacBio 或 Nanopore 數(shù)據(jù)）通常具有較高的錯誤率丁侄，可能會導(dǎo)致組裝的主要序列中存在錯誤和噪聲。直接使用這些序列進行 Hi-C 染色體掛載可能會將這些錯誤傳播到聯(lián)系矩陣中，影響到染色體的準(zhǔn)確性妥箕。

2. 聯(lián)系矩陣質(zhì)量下降： 染色體掛載的關(guān)鍵是基于 Hi-C 數(shù)據(jù)生成染色體之間的聯(lián)系矩陣。如果主要序列存在錯誤更舞，這些錯誤可能會影響到聯(lián)系矩陣的質(zhì)量畦幢，從而降低染色體掛載的準(zhǔn)確性。

3. 組裝斷裂： 主要序列中的錯誤可能會導(dǎo)致組裝斷裂缆蝉，使得 Hi-C 數(shù)據(jù)在某些區(qū)域無法正確匹配宇葱。這會導(dǎo)致染色體掛載失敗或出現(xiàn)錯誤。

4. 重復(fù)序列挑戰(zhàn)： 高復(fù)雜性和重復(fù)性區(qū)域在 de novo 組裝中是一個挑戰(zhàn)刊头，而 CLR 數(shù)據(jù)由于錯誤率較高黍瞧，可能在這些區(qū)域中表現(xiàn)不佳。這會導(dǎo)致在 Hi-C 染色體掛載中難以準(zhǔn)確匹配和定位這些區(qū)域原杂。

因此印颤，為了獲得更準(zhǔn)確的染色體掛載結(jié)果，建議在進行 Hi-C 染色體掛載之前穿肄，對主要序列進行糾錯（polishing）年局，以減少錯誤和噪聲。糾錯可以使用高質(zhì)量的測序數(shù)據(jù)（如 illumina 數(shù)據(jù)）對主要序列進行校正咸产，提高組裝的質(zhì)量和準(zhǔn)確性矢否。這樣，在染色體掛載過程中脑溢，你會使用更可靠的主要序列和聯(lián)系矩陣僵朗，從而獲得更準(zhǔn)確的染色體邊界和結(jié)構(gòu)信息。

另外屑彻，我們可以同時用三代和二代數(shù)據(jù)一起對我們的初始組裝進行糾錯验庙。三代數(shù)據(jù)糾錯主要用于解決長重復(fù)區(qū)域和復(fù)雜結(jié)構(gòu)，提高基因組的連續(xù)性和準(zhǔn)確性酱酬，有助于捕獲染色體級別的信息壶谒；而二代數(shù)據(jù)糾錯主要用于在主要序列中糾正小錯誤，提高基因組的局部準(zhǔn)確性膳沽，從而使基因組序列精準(zhǔn)度更接近于真實的 DNA 序列汗菜。綜合而言，使用三代+二代數(shù)據(jù)進行糾錯是一種綜合優(yōu)化的策略挑社，有助于克服不同測序技術(shù)的優(yōu)缺點陨界，提高基因組的整體準(zhǔn)確性和質(zhì)量。這兩種數(shù)據(jù)的結(jié)合能夠產(chǎn)生更好的組裝和糾錯效果痛阻，為后續(xù)的基因組研究提供更可靠的基因組參考菌瘪。

這里我們我們用了兩種策略：

• 策略1：racon+pilon，即用racon軟件進行三代數(shù)據(jù)糾錯，再用pilon軟件進行二代數(shù)據(jù)糾錯俏扩；
• 策略2：直接使用nextpolish軟件進行三代和二代數(shù)據(jù)同時糾錯糜工。

數(shù)據(jù)準(zhǔn)備

1.primary assembly初始組裝文件，這里用flye組裝的初始組裝做演示录淡。
1.用于的組裝的pacbio CLR三代測序文件捌木。
2.經(jīng)過fastp質(zhì)控過的illumina二代測序文件。

polish策略1

racon+pilon嫉戚，即用racon軟件進行三代數(shù)據(jù)糾錯刨裆，再用pilon軟件進行二代數(shù)據(jù)糾錯。

1# 軟件安裝
racon安裝
軟件主頁：https://github.com/isovic/racon

## 自動安裝
conda install -c bioconda racon

## 手動安裝
git clone --recursive https://github.com/lbcb-sci/racon.git racon
cd racon
mkdir build
cd build
cmake -DCMAKE_BUILD_TYPE=Release ..
make

pilon安裝
軟件主頁：https://github.com/broadinstitute/pilon

## 由于是個java程序直接下載即可使用
wget https://github.com/broadinstitute/pilon/releases/download/v1.24/pilon-1.24.jar 
## java程序統(tǒng)一使用方法：java -jar 后面接java程序
java -jar pilon-1.24.jar --help  # 查看pilon用法

其他軟件依賴：因為涉及到比對彬檀，所以需要安裝三代和二代比對軟件minimap2和bwa
minimap2主頁：https://github.com/lh3/minimap2
bwa主頁：https://github.com/lh3/bwa

## minimap2自動安裝
conda install -c bioconda minimap2
## minimap2手動安裝
wget https://github.com/lh3/minimap2/releases/download/v2.26/minimap2-2.26_x64-linux.tar.bz2
tar -jxvf minimap2-2.26_x64-linux.tar.bz2
./minimap2-2.26_x64-linux/minimap2
---------------------------------------------------------------------------------------------------------------------------------
## bwa手動安裝
git clone https://github.com/lh3/bwa.git
cd bwa
make

2# 運行程序
我這里我打算先用racon進行3輪三代數(shù)據(jù)糾錯帆啃，再用pilon對racon的最后一輪結(jié)果再進行3輪二代數(shù)據(jù)糾錯，加起來相當(dāng)于6輪了窍帝。

• 用racon進行3輪三代數(shù)據(jù)糾錯

## 第1輪糾錯努潘，這里詳細(xì)注釋了，后兩輪就不展開了
/newlustre/home/jfgui/wy/anaconda3/envs/pengfang/bin/minimap2 \  # 用minimap2比對
      -ax map-pb \  # 用-a參數(shù)生成sam文件而非默認(rèn)的paf格式盯桦，-x參數(shù)為比對模式慈俯，這里是pacbio數(shù)據(jù)比對到參考基因組
      -t 24 \ # 線程數(shù)
      ../../purge_dups/flye_purge/assembly.fasta \ #  參考基因組，這里第1輪就是用的flye組裝出的初始組裝拥峦，后面就都是前一輪的糾錯結(jié)果贴膘。
      ../../P01TYD20308306-1_r64030_20201110_065049_1_A02.subreads.fasta \ # 用于糾錯三代測序文件 
      | gzip -c - > minimap1.sam.gz  # 壓縮后輸出為 minimap1.sam.gz
/newlustre/home/jfgui/wy/anaconda3/envs/pengfang/bin/racon \ # 用racon糾錯
      -t 24 \ # 線程數(shù)
      -u \ # 糾錯過和未被就錯的contig都會被輸出，否則只輸出糾錯過的contig
      ../../P01TYD20308306-1_r64030_20201110_065049_1_A02.subreads.fasta \ # 用于糾錯三代測序文件
      minimap1.sam.gz \ # 前面的比對結(jié)果
      ../../purge_dups/flye_purge/assembly.fasta > flye.racon1.fasta # 糾錯前的基因組指向糾錯后的基因組

## 第2輪糾錯
/newlustre/home/jfgui/wy/anaconda3/envs/pengfang/bin/minimap2 -ax map-pb -t 24 flye.racon1.fasta ../../P01TYD20308306-1_r64030_20201110_065049_1_A02.subreads.fasta | gzip -c - > minimap2.sam.gz
/newlustre/home/jfgui/wy/anaconda3/envs/pengfang/bin/racon -t 24 -u ../../P01TYD20308306-1_r64030_20201110_065049_1_A02.subreads.fasta minimap2.sam.gz flye.racon1.fasta > flye.racon2.fasta

## 第3輪糾錯
/newlustre/home/jfgui/wy/anaconda3/envs/pengfang/bin/minimap2 -ax map-pb -t 24 flye.racon2.fasta ../../P01TYD20308306-1_r64030_20201110_065049_1_A02.subreads.fasta | gzip -c - > minimap3.sam.gz
/newlustre/home/jfgui/wy/anaconda3/envs/pengfang/bin/racon -t 24 -u ../../P01TYD20308306-1_r64030_20201110_065049_1_A02.subreads.fasta minimap3.sam.gz flye.racon2.fasta > flye.racon3.fasta

• 再用pilon進行3輪二代數(shù)據(jù)糾錯

## 第1輪糾錯略号，這里詳細(xì)注釋了刑峡，同理后兩輪就不展開了
ln -s ../flye.racon3.fasta  ./racon3.fasta  # 將需要糾錯的基因組鏈接到到當(dāng)前工作路徑
/newlustre/home/jfgui/Wangtao/software/bwa/bwa index racon3.fasta    # 對基因組構(gòu)建index以便后面的比對
/newlustre/home/jfgui/Wangtao/software/bwa/bwa mem \ # 用bwa的mem比對模式對二代數(shù)據(jù)進行比對
    -t 24 \ #線程數(shù)
    racon3.fasta \ # 對比用的參考基因組，即前一輪的糾錯結(jié)果
    ../../../../../XX_survey/fastp.out/liver.clean.CLEAN.R1.fastq.gz ../../../../../XX_survey/fastp.out/liver.clean.CLEAN.R2.fastq.gz \ # 用于糾錯的二代數(shù)據(jù)
    | /newlustre/home/jfgui/wy/anaconda3/envs/pengfang/bin/samtools sort \ # 用samtools的sort命令對比對結(jié)果進行排序玄柠，用管道符連接 
           -@ 24 \ # samtools的線程數(shù)
           -  \ # 表示管道符前面命令的輸出作為管道符后面命令的輸入突梦，這里也就是指未排序的bam文件
           -o bwamem1.bam # 輸出比對結(jié)果
/newlustre/home/jfgui/wy/anaconda3/envs/pengfang/bin/samtools index bwamem1.bam # 用samtools對比對文件進行index
java -Xmx128G -jar /newlustre/home/jfgui/Wangtao/software/pilon-1.24.jar \  # 用pilon糾錯，-Xmx128G指的是用128G內(nèi)存
    --genome racon3.fasta \  # 前一輪糾錯的基因組
    --frags bwamem1.bam \   # 比對結(jié)果文件
    --changes \ # 糾錯前后變化記錄文件
    --diploid \ # 二倍體羽利，不設(shè)置會默認(rèn)單倍體進而識別單倍型
    --threads 24 \   # 線程數(shù)
    --outdir ./pilon1.out \  # 輸出目錄 
    --output racon3.pilon1   # 輸出前綴

## 第2輪糾錯
ln -s ./pilon1.out/racon3.pilon1.fasta  ./racon3.pilon1.fasta
/newlustre/home/jfgui/Wangtao/software/bwa/bwa index racon3.pilon1.fasta
/newlustre/home/jfgui/Wangtao/software/bwa/bwa mem -t 24 racon3.pilon1.fasta ../../../../../XX_survey/fastp.out/liver.clean.CLEAN.R1.fastq.gz ../../../../../XX_survey/fastp.out/liver.clean.CLEAN.R2.fastq.gz | /newlustre/home/jfgui/wy/anaconda3/envs/pengfang/bin/samtools sort -@ 24 - -o bwamem2.bam
/newlustre/home/jfgui/wy/anaconda3/envs/pengfang/bin/samtools index bwamem2.bam
java -Xmx128G -jar /newlustre/home/jfgui/Wangtao/software/pilon-1.24.jar --genome racon3.pilon1.fasta --frags bwamem2.bam --changes  --diploid --threads 24 --outdir ./pilon2.out --output racon3.pilon2

## 第3輪糾錯
ln -s ./pilon2.out/racon3.pilon2.fasta  ./racon3.pilon2.fasta
/newlustre/home/jfgui/Wangtao/software/bwa/bwa index racon3.pilon2.fasta
/newlustre/home/jfgui/Wangtao/software/bwa/bwa mem -t 24 racon3.pilon2.fasta ../../../../../XX_survey/fastp.out/liver.clean.CLEAN.R1.fastq.gz ../../../../../XX_survey/fastp.out/liver.clean.CLEAN.R2.fastq.gz | /newlustre/home/jfgui/wy/anaconda3/envs/pengfang/bin/samtools sort -@ 24 - -o bwamem3.bam
/newlustre/home/jfgui/wy/anaconda3/envs/pengfang/bin/samtools index bwamem3.bam
java -Xmx128G -jar /newlustre/home/jfgui/Wangtao/software/pilon-1.24.jar --genome racon3.pilon2.fasta --frags bwamem3.bam --changes  --diploid --threads 24 --outdir ./pilon3.out --output racon3.pilon3

3# 結(jié)果查看

3輪racon結(jié)果

$assembly-stats flye.racon1.fasta flye.racon2.fasta flye.racon3.fasta
stats for flye.racon1.fasta
sum = 774215185, n = 680, ave = 1138551.74, largest = 25240035
N50 = 9242029, n = 27
N60 = 7581968, n = 36
N70 = 6421899, n = 47
N80 = 5066544, n = 61
N90 = 2268818, n = 81
N100 = 517, n = 680
N_count = 0
Gaps = 0
-------------------------------------------------------------------------------
stats for flye.racon2.fasta
sum = 775085635, n = 680, ave = 1139831.82, largest = 25259389
N50 = 9253497, n = 27
N60 = 7586889, n = 36
N70 = 6427024, n = 47
N80 = 5073214, n = 61
N90 = 2269709, n = 81
N100 = 517, n = 680
N_count = 0
Gaps = 0
-------------------------------------------------------------------------------
stats for flye.racon3.fasta
sum = 775642180, n = 680, ave = 1140650.26, largest = 25269729
N50 = 9260756, n = 27
N60 = 7588618, n = 36
N70 = 6429082, n = 47
N80 = 5078123, n = 61
N90 = 2270025, n = 81
N100 = 517, n = 680
N_count = 0
Gaps = 0

pilon跑的非常慢！不知道是不是哪里參數(shù)設(shè)置得不對匾浪，3輪pilon糾錯跑了大概2.5天皇帮，結(jié)果如下：

pilon輸出結(jié)果

最終糾錯結(jié)果分別在pilon1.out，pilon2.out蛋辈，pilon3.out中属拾，每個目錄里包含了糾錯后的基因組fasta文件及糾錯過程changes文件。
assembly-stats查看，基因組大小依次變小渐白，準(zhǔn)確性應(yīng)該是上升了尊浓，等待后續(xù)組裝質(zhì)量評估。

$assembly-stats pilon1.out/racon3.pilon1.fasta pilon2.out/racon3.pilon2.fasta pilon3.out/racon3.pilon3.fasta
stats for pilon1.out/racon3.pilon1.fasta
sum = 774270266, n = 680, ave = 1138632.74, largest = 25232310
N50 = 9243114, n = 27
N60 = 7577994, n = 36
N70 = 6420177, n = 47
N80 = 5068385, n = 61
N90 = 2266107, n = 81
N100 = 517, n = 680
N_count = 0
Gaps = 0
-------------------------------------------------------------------------------
stats for pilon2.out/racon3.pilon2.fasta
sum = 774210928, n = 680, ave = 1138545.48, largest = 25230882
N50 = 9242790, n = 27
N60 = 7578144, n = 36
N70 = 6419986, n = 47
N80 = 5068374, n = 61
N90 = 2266189, n = 81
N100 = 517, n = 680
N_count = 0
Gaps = 0
-------------------------------------------------------------------------------
stats for pilon3.out/racon3.pilon3.fasta
sum = 774169945, n = 680, ave = 1138485.21, largest = 25229904
N50 = 9242646, n = 27
N60 = 7577699, n = 36
N70 = 6420023, n = 47
N80 = 5068254, n = 61
N90 = 2266164, n = 81
N100 = 517, n = 680
N_count = 0
Gaps = 0

polish策略2

直接使用nextpolish軟件進行三代和二代數(shù)據(jù)同時糾錯纯衍。
1# 軟件安裝
nextpolish安裝
軟件主頁：https://github.com/Nextomics/NextPolish

## 手動安裝
pip install paralleltask
wget https://github.com/Nextomics/NextPolish/releases/download/v1.4.1/NextPolish.tgz
tar -vxzf NextPolish.tgz
cd NextPolish
make

## 自動安裝
conda install -c bioconda nextpolish

2# 運行程序
nextPolish的使用跟前面兩個軟件有些不同眠砾，是要先填寫一個.cfg配置文件，然后根據(jù)這個配置文件軟件來決定怎么開展polish托酸。軟件官網(wǎng)提供了4種配置文件案例，分別為只用二代短reads糾柒巫、只用三代長reads糾励堡、同時用短reads和長reads糾以及用短reads和hifi reads糾。我們這次選擇第三種配置方式堡掏。

• 首先準(zhǔn)備短reads路徑文件

ls ../../../../XX_survey/fastp.out/liver.clean.CLEAN.R1.fastq.gz  ../../../../XX_survey/fastp.out/liver.clean.CLEAN.R2.fastq.gz > nextpolish.sgs.fofn

• 同時準(zhǔn)備長reads路徑文件

ls ../../P01TYD20308306-1_r64030_20201110_065049_1_A02.subreads.fasta > nextpolish.lgs.fofn

• 其次配置.cfg文件

$ vi nextpolish.run.cfg
------------------------------------------------------------------------------------------------
[General]
job_type = local     # 運行方式本地還是集群（local, sge, pbs）
job_prefix = nextPolish
task = best     # 運行策略应结，2輪三代糾錯+4輪二代糾錯，如果后面只有三代或者二代數(shù)據(jù)泉唁，則只運行三代或者二代糾錯
rewrite = yes
rerun = 3
parallel_jobs = 6    # 平行任務(wù)數(shù)
multithread_jobs = 4    # 每個任務(wù)所給的線程數(shù)
genome = ../../purge_dups/flye_purge/assembly.fasta    # 初始基因組
genome_size = auto    # 自動計算上面給的基因組大小
workdir = ./nextpolish_rundir    # 建工作路徑
polish_options = -p {multithread_jobs}     # 按照上面配置的平行任務(wù)的策略跑

[sgs_option]
sgs_fofn = ./nextpolish.sgs.fofn    # 前面生成的二代測序數(shù)據(jù)路徑文件
sgs_options = -max_depth 100 -bwa     # bwa比對鹅龄，用最大100*的數(shù)據(jù)去做糾錯，也可以用全部的

[lgs_option]
lgs_fofn = ./nextpolish.lgs.fofn     # 前面生成的二三代測序數(shù)據(jù)路徑文件
lgs_options = -min_read_len 1k -max_depth 100     #minimap2比對亭畜，同樣用最大100*的數(shù)據(jù)去做糾錯
lgs_minimap2_options = -x map-pb      #minimap2比對方式扮休，pacbio數(shù)據(jù)比對到基因組上。

• 最后nextPolish運行.cfg即可

/newlustre/home/jfgui/Wangtao/software/NextPolish/nextPolish  ./nextpolish.run.cfg

3# 結(jié)果查看

nextpolish輸出結(jié)果

最后我們查看糾錯前后結(jié)果對比汰现，基因組大小糾錯后少了0.5M挂谍，準(zhǔn)確性肯定提高了只不過這個方法看不出來??

$assembly-stats  nextpolish_rundir/00.lgs_polish/input.genome.fasta  nextpolish_rundir/genome.nextpolish.fasta
stats for nextpolish_rundir/00.lgs_polish/input.genome.fasta # 糾錯前
sum = 771989320, n = 680, ave = 1135278.41, largest = 25186124
N50 = 9211426, n = 27
N60 = 7566198, n = 36
N70 = 6405554, n = 47
N80 = 5050486, n = 61
N90 = 2266136, n = 81
N100 = 493, n = 680
N_count = 0
Gaps = 0
-------------------------------------------------------------------------------
stats for nextpolish_rundir/genome.nextpolish.fasta  # 糾錯后
sum = 771468932, n = 680, ave = 1134513.14, largest = 25167817 
N50 = 9204720, n = 27
N60 = 7567313, n = 36
N70 = 6408529, n = 47
N80 = 5048648, n = 61
N90 = 2265774, n = 81
N100 = 493, n = 680
N_count = 0
Gaps = 0

其他相關(guān)好文推薦：racon+pilon進行三代+二代數(shù)據(jù)糾錯 - 簡書 (jianshu.com)