结构变异（SV）鉴定软件——cuteSV(一)

0.cuteSV简介

cuteSV是一款利用三代Longreads比对鉴定结构变异（SV）的软件，于2020发布在Genome Biol上。

github：tjiangHIT/cuteSV: Long read based human genomic structural variation detection with cuteSV (github.com)

优点是快速准确且易安装使用。个人感觉是几款基于reads比对检测SV的软件中最好用的一款。

此外，cuteSV也支持基于全基因组比对call SV，可见 ：结构变异（SV）鉴定软件——cuteSV(二)-CSDN博客

1.cuteSV安装

安装很简单，可以直接用pip安装

pip install cuteSV

如果有什么问题，pip装不上也可以用conda

conda install -c bioconda cutesv

2.reads比对

作为鉴定SV的软件，第一步需要long reads比对到参考基因组上。cuteSV官方示例的比对软件是minimap2。

#Hi-Fi(CCS)数据
minimap2 --paf-no-hit -ax map-hifi -t 30 ref.fasta pacbio-ccs.fq.gz |samtools view -bS - |samtools sort -@ 30 -o minimap2.bam

#ONT数据
minimap2 --paf-no-hit -ax map-ont -t 30 ref.fa ont.fq.gz | samtools view -bS - |samtools sort -@ 30 -o minimap2.bam

#CLR数据
minimap2 --paf-no-hit -ax map-pb -t 30 ref.fa pacbio.fq.gz | samtools view -bS - |samtools sort -@ 30 -o minimap2.bam

 之后对生成的bam文件建立索引

samtools index -@ 30 minimap2.bam

3.运行cuteSV

首先是一些重要的参数：

-q 最低MIN_MAPQ得分，[默认20]

-l 最小SV长度，[默认30]

-L 最大SV长度, [默认 100000]

-s 最小reads覆盖度，[默认10]

--genotype 区分单倍型

此外还有一些高级的参数，cuteSV已经给出了针对不同的reads类型给出了推荐参数：

> For PacBio CLR data:
	--max_cluster_bias_INS		100
	--diff_ratio_merging_INS	0.3
	--max_cluster_bias_DEL	200
	--diff_ratio_merging_DEL	0.5

> For PacBio CCS(HIFI) data:
	--max_cluster_bias_INS		1000
	--diff_ratio_merging_INS	0.9
	--max_cluster_bias_DEL	1000
	--diff_ratio_merging_DEL	0.5

> For ONT data:
	--max_cluster_bias_INS		100
	--diff_ratio_merging_INS	0.3
	--max_cluster_bias_DEL	100
	--diff_ratio_merging_DEL	0.3
> For force calling:
	--min_mapq 			10

以下为使用上述推荐高级参数，-l 设为 50 -s 设为5，其余均为默认参数的示例；其中$accession指的是样品名，$work_dir 是输出目录。

CLR数据：

cuteSV --genotype -l 50 -s 5 -S $accession --max_cluster_bias_INS 100 --diff_ratio_merging_INS 0.3 --max_cluster_bias_DEL 200 --diff_ratio_merging_DEL 0.5 $accession.minimap2.bam ref.fasta accession.cuteSV.vcf $work_dir

CCS（HIFI）数据：

cuteSV --genotype -l 50 -s 5 -S $accession --max_cluster_bias_INS 1000 --diff_ratio_merging_INS 0.9 --max_cluster_bias_DEL 1000 --diff_ratio_merging_DEL 0.5 $accession.minimap2.bam ref.fasta accession.cuteSV.vcf $work_dir

ONT数据：

cuteSV --genotype -l 50 -s 5 -S $accession --max_cluster_bias_INS 100 --diff_ratio_merging_INS 0.3 --max_cluster_bias_DEL 100 --diff_ratio_merging_DEL 0.3 $accession.minimap2.bam ref.fasta accession.cuteSV.vcf $work_dir

4.提供一个pipeline

在这里提供一个针对HIFI reads 比对加call SV的pipeline。使用方法为：

perl cuteSV_pipline.pl -r reference.fasta -q quary.fasta -a AJ -w work_dir

-r 参考基因组 fasta文件
-q HIFI reads 文件（.fa 或者 .fq）
-a vcf样品前缀
-w 输出目录地址

具体代码如下：

use strict;
use warnings;
use Getopt::Long;
my ($ref,$help,$accession,$quary);

GetOptions(
        "ref=s" =>\$ref,
        "accession=s" =>\$accession,
		"quary=s" =>\$quary,
        "help+" =>\$help,
		"work=s" =>\$wd
);




system("mkdir -p $wd");
system("cd $wd");
#minimap2
$cmdString ="minimap2 --paf-no-hit -ax map-hifi -t 30 $ref $quary |samtools view -bS - |samtools sort -@ 30 -o $wd/$accession.minimap2.bam -";
print STDOUT (localtime) . ": CMD: $cmdString\n";
 system("$cmdString") == 0 or die "failed to execute: $cmdString\n";

#index
$cmdString ="samtools index -@ 30 $wd/$accession.minimap2.bam";
print STDOUT (localtime) . ": CMD: $cmdString\n";
system("$cmdString") == 0 or die "failed to execute: $cmdString\n";

#cuteSV
$cmdString ="cuteSV --genotype -l 50 -s 5 -S $accession --max_cluster_bias_INS 1000 --diff_ratio_merging_INS 0.9 --max_cluster_bias_DEL 1000 --diff_ratio_merging_DEL 0.5 $wd/$accession.minimap2.bam $ref $wd/$accession.cuteSV.vcf $wd";
print STDOUT (localtime) . ": CMD: $cmdString\n";
system("$cmdString") == 0 or die "failed to execute: $cmdString\n";