寫在前面
- 以下內(nèi)容均來自我在菲沙基因（Frasergen）暑期生信培訓(xùn)班上記錄的課堂筆記
1.Hi-C 輔助組裝（PGA）技術(shù)原理
染色質(zhì)在細(xì)胞核內(nèi)分布的并不是隨機分布的，而是不同染色體占據(jù)不同的空間。
-
Hi-C實驗原理圖
-
基因組互作衰減
染色體內(nèi)互作強度較強裂七，但也隨著空間距離的增大互作強度在衰減
-
Hi-C互作三大規(guī)律
- 1.染色體內(nèi)互作富集（Intrachromosomal interaction enrichment）
- 2.互作隨距離衰減（distance-dependent interaction decay）
- 3.局部互作平滑（local interaction smoothness）
2.Hi-C技術(shù)輔助基因組組裝目的
- 將contigs組裝到假染色體層面，實現(xiàn)基因組組裝到染色體層面
3.Hi-C技術(shù)主流輔助組裝軟件
-
首選軟件：Lachesis
-
輔助組裝步驟
1.Cluster(聚類)龙誊。因為染色質(zhì)內(nèi)的互 作強度要高于染色質(zhì)間的互作強度，所 以先對contig/scaffold進行聚類成染色體群喷楣。
2.Order(排序)载迄。確定每個染色體群中contig/scaffold的順序
3.Orient(定向)。確定每個 contig/scaffold的方向三個步驟按照互作強度依次
-
Cluster聚類原理
-
Order排序原理
-
Orient定向原理
4.Hi-C PGA代碼實操(PGA即輔助組裝)
4.1 軟件安裝與編譯
git clon https://github.com/shendurelab/LACHESIS.git
#此外還需要使用conda安裝boost與samtools
conda install -c conda-forge boost-cpp -y
conda install -c bioconda samtools -y
conda install -c bioconda bwa -y
#都安裝完成后抡蛙，開始編譯
cd LACHESIS/
make
#后續(xù)用到的腳本都在/lachesis/src/bin目錄下
4.2 數(shù)據(jù)比對
- 在基因組草圖目錄中建立bwa索引
bwa index draft_genome.fasta
- 采用bwa-mem模式，將Hi-C reads比對到基因組草圖魂迄，生成sam文件
bwa mem -t 10 ./draft_genome.fasta ./R1.fastq.gz ./R2.fastq.gz > thaliana.sam
4.3 數(shù)據(jù)過濾
- 使用lachesis首先過濾掉flagstat大于2048的序列(flagstat > 2048表示去除補充匹配的reads見下圖)
perl /lachesis/PreprocessSAMs-rmsubalig.pl thaliana.sam thaliana.II.sam
#PreprocessSAMs-rmsubalig.pl腳本接輸入和輸出文件的名字粗截，運行完會得到一個過濾后的sam文件：后綴.II.sam
- 將PreprocessSAMs.pl程序復(fù)制到當(dāng)前目錄，并加以修改
cp /lachesis/src/bin/PreprocessSAMs.pl ./
vim PreprocessSAMs.pl
#小tips：vim下使用set nu顯示行號捣炬，輸入數(shù)字直接跳轉(zhuǎn)
#根據(jù)Hi-C實驗過程中使用的限制性內(nèi)切酶熊昌，更改酶切位點序列 (擬南芥使用的是HindIII，對應(yīng)的酶切位點序列是AAGCTT湿酸，$RE_site對應(yīng)的就是酶切序列)
- 將PreprocessSAMs.sh程序復(fù)制到當(dāng)前目錄婿屹，并加以修改
cp /lachesis/src/bin/PreprocessSAMs.sh ./
vim PreprocessSAMs.sh
#以下幾處需修改：
#DIR
#在第46行，修改為目前所在的目錄路徑推溃，可用pwd命令顯示當(dāng)前路徑
#SAMS
#在第48行昂利，修改為上一步生成的sam文件名字 (例如thaliana.II.sam)
#ASSEMBLY
#在第49行，修改為參考基因組的路徑(例如../01.ref/arabidopsis_thaliana_draft_genome.fasta)
#注意：ASSEMBLY一定要用相對路徑，否則會報錯PreprocessSAMs.pl: Can't find draft assembly fil
- 然后過濾出距離酶切位點500bp以內(nèi)且唯一比對的reads
sh PreprocessSAMs.sh
# 注意：這一步運行完需要按回車鍵結(jié)束
cd ..
- 結(jié)果文件
4.4 Lachesis組裝
- 建立一個lachesis分析文件夾并進入蜂奸，例如：
mkdir 04.lachesis
cd 04.lachesis
- 新建一個過濾文件的目錄犁苏，同時建立過濾文件的軟鏈接
mkdir bam_file
cd bam_file
ln -s /data/alnbam/*.II.REduced.paired_only.bam
cd ..
- 復(fù)制一個lachesis的配置文件到當(dāng)前目錄，并修改
cp lachesis/src/bin/INIs/test_case.ini ./
vim test_case.ini
-
以下為輸入輸出文件位置參數(shù)扩所，本例需修改：
OUTPUT_DIR：第35行围详，建議修改為out/NNJ_90_2_3_120_10，用文件夾名來記錄重要參數(shù)值
DRAFT_ASSEMBLY_FASTA：第45行祖屏，修改為參考基因組的路徑助赞，建議設(shè)為絕對路徑/local_data1/draft_genome.fasta
SAM_DIR：第49行，修改為比對文件所在的路徑袁勺，即bam_file
SAM_FILES：第54行雹食，修改為過濾得到的bam文件的名字，本例中擬南芥有兩組數(shù)據(jù)：samplex.II.REduced.paired_only.bam sampley.II.REduced.paired_only.bam
RE_SITE_SEQ：第58行魁兼，因為本例中是使用的是HindIII酶婉徘，其對應(yīng)的序列為AAGCTT，因此酶切序列設(shè)為AAGCTT
修改test_case.ini -
配置文件test_case.ini的重要參數(shù)咐汞。在實際項目中會有不同
USE_REFERENCE：0代表不使用標(biāo)準(zhǔn)基因組進行評估盖呼；1表示用標(biāo)準(zhǔn)基因組進行評估（本例中不用修改）
REF_ASSEMBLY_FASTA：后面是標(biāo)準(zhǔn)參考基因組，如：REF_ASSEMBLY_FASTA = /local_data1/ninanjie.fa
BLAST_FILE_HEAD：后面加標(biāo)準(zhǔn)基因組和參考基因組的blastn比對結(jié)果化撕，這里的格式固定几晤，只需輸入前綴BLAST_FILE_HEAD = /local_data1/ler1_vs_tari10CLUSTER_N：第129行，后面接物種所對應(yīng)的染色體數(shù)目植阴，本例樣本為擬南芥蟹瘾，因此設(shè)為5
CLUSTER_MIN_RE_SITES：第137行，聚類分析中contig/scaffold序列中的最小酶切位點數(shù)掠手，建議設(shè)為90
CLUSTER_MAX_LINK_DENSITY：第140行憾朴，聚類分析中contig/scaffold序列中的最大Link深度，建議設(shè)為2
CLUSTER_NONINFORMATIVE_RATIO：第144行喷鸽，聚類回插中contig/scaffold序列與目標(biāo)cluster互作數(shù)同其他cluster互作數(shù)比例众雷，建議設(shè)為3
ORDER_MIN_N_RES_IN_TRUNK：第153行，排序分析中TRUNK內(nèi)contig/scaffold序列中的最小酶切位點數(shù)做祝，建議設(shè)為120
ORDER_MIN_N_RES_IN_SHREDS：第155行砾省，排序分析中回插contig/scaffold序列的最小酶切位點數(shù)，建議設(shè)為10
最后新建的兩個目錄作為組裝基因組和標(biāo)準(zhǔn)基因組比對圖片的存放地址
mkdir ~/bin #會在這個目錄下生成畫圖腳本
mkdir ~/public_html #在這個目錄下生成圖片
ulimit -s unlimited #表示不限制線程混槐，如果沒有設(shè)置則會報錯
-
運行Lachesis
/gnu_modulefile/lachesis test_case.ini
- 結(jié)果文件
txt文件是Lachesis組裝的基因和已發(fā)表的標(biāo)準(zhǔn)基因組的比較結(jié)果信息编兄，Lachesis的組裝信息在out文件夾中進入我們設(shè)置的out/NNJ_90_2_3_120_10
main_results/中的clusters.txt和clusters.by_name.txt展示了lachesis聚類的結(jié)果，一行代表一個group（相當(dāng)于一條染色體）group*.ordering則展示的是每一個group中排序和方向的結(jié)果
cached_data/中則是各個group互作矩陣文件声登，以及在lachesis排序過程中的trunk信息狠鸳。這個文件夾存放的是程序運行的中間臨時文件揣苏，僅知道即可
REPORT.txt則是對Lachesis的組裝結(jié)果做了一個總體的展示.REPORT.txt第一部分記錄Lachesis組裝全部參數(shù)；REPORT.txt第二部分記錄Cluster碰煌，Order以及最后Trunk的信息：
CLUSTERING METRICS部分：如果給了近緣物種信息舒岸，會將Cluster與近緣物種進行對比，統(tǒng)計沒有比對上的Contig數(shù)以及在染色體上位置不同的Contig數(shù)芦圾；如果項目沒有提供近緣物種信息蛾派，則僅統(tǒng)計Cluster信息
ORDERING METRICS部分：記錄了orderings中contigs數(shù)目以及堿基數(shù)，trunks中的contigs數(shù)目以及堿基數(shù)个少。其中orderings中contigs長度與所有序列長度的比值即掛載率洪乍，如下圖中紅色框內(nèi)的
- 通過讀取main_results文件夾中的聚類、排序與定向信息夜焦，可以使用lachesis中CreateScaffoldedFasta.pl腳本生成最終的組裝結(jié)果fasta文件
perl /src/bin/CreateScaffoldedFasta.pl /local_data1/draft_genome.fasta /local_data1/04.lachesis/out/NNJ_90_2_3_120_10
#輸入兩個參數(shù)壳澳，第一個是參考基因組路徑，第二個是上一步lachesis組裝時的輸出文件夾路徑
- 生成的最終組裝文件：
/local_data1//04.lachesis/out/NNJ_90_2_3_120_10/Lachesis_assembly.fasta
