这个软件好像都没有什么人用哎~
但是我今天用上了，为什么会用上呢，还不是在KB的组装那个三倍体的基因组，为了区分测序回来的那么那么多的reads都是属于谁家的！！！！
对于这个软件，需要一些依赖，不过基本上这些依赖我们早就用过了，也就不存在依赖不依赖的了吧~
需要的依赖环境
浅浅列举一下吧：做到这一步，这些软件还没有的话，有点说不过去吧？
1. MUMmer Version:4.0.0beta2 ([https://github.com/mummer4/mummer](https://github.com/mummer4/mummer))
2. Bismark Version:v0.22.1 ([https://github.com/FelixKrueger/Bismark](https://github.com/FelixKrueger/Bismark))
3. minimap2 Version:2.17 ([https://github.com/lh3/minimap2](https://github.com/lh3/minimap2))
4. HISAT2 Version:2.1.0 ([http://ccb.jhu.edu/software/hisat2/manual.shtml](http://ccb.jhu.edu/software/hisat2/manual.shtml))
5. pysam Version:0.15.2 ([https://pypi.org/project/pysam/](https://pypi.org/project/pysam/))
6. SAMtools Version:1.6 ([http://samtools.sourceforge.net/](http://samtools.sourceforge.net/))
Installation:
git clone https://github.com/jia-wu-feng/PP2PG.git
cd PP2PG
chmod +X tool/phasing.py
这是github上的一张流程介绍
Fig.1 The phasing pipeline based on two parental genomes (PP2PG). It mainly consists of four steps: (i) Mapping, (ii) SNP Calling, (iii) Phasing and (iv) Classification.
也就是包括四个步骤：（i）映射，（ii）SNP calling，（iii）评分和（iv）分类。
(i)mapping(四选一)
1.Iso-Seq /PicBio读取的类型使用 minimap2 软件进行对齐：
minimap2 -ax splice -uf --secondary=no -C5 -O6,24 -B4 --MD ref query > results.sam
2.利用hisat2软件对RNA-Seq配对双末端读段的类型进行比对：
hisat2 --dta -k 1 -x genome_index -1 query_1 -2 query_2 -S out.sam
3.BS-Seq读取的类型使用Bismark软件进行对齐：
bismark --bowtie2 -genome genome_index -1 query_1 -2 query_2
4.Ribo-seq读取的类型使用Tophat2（Bowtie2）软件进行对齐：
tophat2 -N 2 -I 50000 -G gff_file/gtf_file -o prefix genome_index query
(ii) SNP Calling
nucmer --mum --maxgap=500 --mincluster=100 --prefix=g1_g2 g1_genome.fa g1_genome.fa
delta-filter -1 -q -r g1_g2.delta > g1_g2.filter.delta
show-snps -C -H -I -T -r -l g1_g2.filter.delta > g1_g2.snp
(iii) Scoring and Classification:
phasing.py --pt=<phase_type> --g1=<genome1_alignment> --g2=<genome2_alignment>
--snp=<snpfile> --gop1=<genome1_output_prefix> --gop2=<genome2_output_prefix>
4.Output file
Two genome-based alignment results file, containing 2 folders. One is the phasing results of g1 genome alignment, which includes: g1-synteny reads, g2-synteny reads, unknown reads, and reads that only mapped to the g1 genome (g1-only reads); the other is the phasing results of g2 genome alignment, which includes: g1-synteny reads, g2-synteny reads, unknown reads, and reads that only mapped to the g2 genome (g2-only reads).
两个基于基因组的比对结果文件，包含 2 个文件夹。 一个是g1基因组比对的阶段结果，其中包括：g1-synteny读段，g2-synteny读段，未知读段和仅映射到g1基因组的读段（仅g1读段）;另一种是G2基因组比对的阶段结果，其中包括：G1-Synteny读段，G2-Synteny读段，未知读段和仅映射到G2基因组的读段（仅G2读段）。
g1_genome/g1reads.bam --- g1-synteny reads mapped to g1 genome
g1_genome/g2reads.bam --- g2-synteny reads mapped to g1 genome
g1_genome/unknowreads.bam --- unknown reads mapped to g1 genome
g1_genome/g1onlyreads.bam --- reads that only mapped to the g1 genome (g1-only reads)
g2_genome/g1reads.bam --- g1-synteny reads mapped to g2 genome
g2_genome/g2reads.bam --- g2-synteny reads mapped to g2 genome
g2_genome/unknowreads.bam --- unknown reads mapped to g2 genome
g2_genome/g2onlyreads.bam --- reads that only mapped to the g2 genome (g2-only reads)
g1_g2_Phasing_Report.txt --- Final Phasing Reads Report
g1_g2_Phasing_Report.txt ---最终阶段阅读报告
Example
==========================
Final Phasing Reads Report
==========================
Phasing Data Type: Iso-Seq (PacBio Isoform Sequence).
Note: When calculating the Separation Rate, the reads which are unmapped in two parental genomes are discarded.
Type of reads: Read Counts
MH63-synteny reads: 129
ZS97-synteny reads: 129
Unknown reads: 47
MH63-only reads: 2
ZS97-only reads: 3
Separation Rate: 83.23%
还是看看自己是怎么用的吧~~~
# 先划分参考的序列，
minimap2 -a -t 40 -x asm10 -m 50 --secondary=no A.contig.mmi m64189e_221213_192604.ccs.fasta >A.sam
minimap2 -a -t 40 -x asm10 -m 50 --secondary=no BB.contig.mmi m64189e_221213_192604.ccs.fasta >BB.sam
#sam文件实在是太占内存了，还是转为bam吧
samtools sort -@ 20 -O bam -o A.sorted.bam A.sam
samtools sort -@ 20 -O bam -o BB.sorted.bam BB.sam
#我就有比对质量高的reads的id
awk '$5 > 55 {print $1}' A.sorted.bam | sort -u > A.sorted.id
#把A和BB的sorted的id先取个交集
cat A.sorted.id BB.sorted.id | sort | uniq -c | awk '$1 == 2{print $2}' > A_BB.JJ.id
#再来取个差集
sort A.sorted.id A_BB.JJ.id A_BB.JJ.id > A.uniq.id
##all.reads.id - A_and_BB.uniq.id(取差集)
sort all.reads.id A_and_BB.uniq.id A_and_BB.uniq.id | uniq -u > unkown.reads.id
#提取unkown的reads序列
seqtk subseq ../ccs.fasta unkown.reads.id > unkown.reads.fa
上面的和这个软件没有什么关系，我就单纯想写，下面才是这个软件怎么用~~~
#step1.Mapping
minimap2 -ax splice -uf --secondary=no -C5 -O6,24 -B4 --MD -t 60 ../A.contig.fa unkown.reads.fa |samtools sort -@ 60 -O bam -o A.sorted.bam -
minimap2 -ax splice -uf --secondary=no -C5 -O6,24 -B4 --MD -t 60 ../BB.contig.fa unkown.reads.fa |samtools sort -@ 60 -O bam -o B.sorted.bam -
#Step2.SNP calling
nucmer --mum --maxgap=500 --mincluster=100 --prefix=A_BB -t 60 ../A.contig.fa ./BB.contig.fa
delta-filter -1 -q -r A_BB.delta > A_BB.filter.delta
show-snps -C -H -I -T -r -l A_BB.filter.delta > A_BB.snp
#step3.Scoring and Classification
/share/home/off/biosoft/PP2PG/tool/phasing.py --pt=Isoseq --g1=A.sorted.bam --g2=B.sorted.bam --snp=A_BB.snp --gop1=A.reads --gop2=B.reads --st 60
我也在等结果！！！
