Control-FREEC的使用

评论910

Control-FREEC是一款广受好评的WGS、WES数据CNV检测软件。

安装

当前最新版本是v11.6（29 May 2020）。

下载源码编译安装

  1. wget https://github.com/BoevaLab/FREEC/archive/refs/tags/v11.6.tar.gz
  2. tar -zxvf v11.6.tar.gz
  3. cd FREEC-11.6/src
  4. make

获得编译后的freec执行程序。

如果需要统计GC比例，还需要安装gccount

  1. wget http://bioinfo-out.curie.fr/projects/freec/src/gccount.tar.gz
  2. tar -zxvf gccount.tar.gz -C gccount

gccount已预编译好，解压后可直接使用。

另外，control-freec需要环境中包含（或配置文件中指定路径）R、samtools、bedtools、sambamba等。其中samtools和sambamba用于处理bam文件，bedtools用于生成minipileup格式文件。

如果需要结果中包含比对信息，需下载mappability track，官方提供hg19、hg38、mm9、mm10等基因组的track。

如果需要检测等位染色体状态，需下载dbsnp文件，官方提供了hg19、hg18、mm10等，不知道和ncbi提供的dbsnp格式是否存在差异。

使用

Control-Freec的使用命令非常简单

  1. freec -conf conf.txt

重要的是配置文件，软件源码库中包含了WES和WGS等两个例子。但是，例子是的内容是相当不全的，更多的配置标签需要看官网的教程。

以下是WGS例子，我增加了输出路径等及删除了部分用不上的参数。

  1. [general]
  2. outputDir = /path/to/output
  3. chrLenFile = b37.freec.len
  4. ploidy = 2
  5.  
  6. breakPointThreshold = .8
  7.  
  8. #coefficientOfVariation = 0.01
  9. window = 50000
  10. #step=10000
  11.  
  12. chrFiles = /path/to/b37/sep_chrom
  13. GCcontentProfile = /path/to/b37/b37.freec.50kb.gc.cnp
  14.  
  15. maxThreads = 8
  16.  
  17. #readCountThreshold=10
  18. #numberOfProcesses = 4
  19. #outputDir = test
  20. #contaminationAdjustment = TRUE
  21. #contamination = 0.4
  22. #minMappabilityPerWindow = 0.95
  23.  
  24. #breakPointType = 4
  25. #forceGCcontentNormalization = 0
  26. sex=XY
  27.  
  28. #BedGraphOutput=TRUE
  29.  
  30. sambamba = sambamba
  31. SambambaThreads = 8
  32.  
  33. [sample]
  34. mateFile = test.bam
  35. #mateCopyNumberFile = test/sample.cpn
  36. inputFormat = BAM
  37. mateOrientation = 0
  38. ##use "mateOrientation=0" for sorted .SAM and .BAM
  39.  
  40. [control]
  41. #mateFile = /path/control.pileup.gz
  42. #mateCopyNumberFile = path/control.cpn
  43. #inputFormat = pileup
  44.  
  45. #mateOrientation = RF
  46.  
  47. #[BAF]
  48. ##use the following options to calculate B allele frequency profiles and genotype status. This option can only be used if "inputFormat=pileup"
  49.  
  50. #SNPfile = /bioinfo/users/vboeva/Desktop/annotations/hg19_snp131.SingleDiNucl.1based.txt
  51. #minimalCoveragePerPosition = 5
  52.  
  53. ##use "minimalQualityPerPosition" and "shiftInQuality" to consider only high quality position in calculation of allelic frequencies (this option significantly slows down reading of .pileup)
  54.  
  55. #minimalQualityPerPosition = 5
  56. #shiftInQuality = 33
  57.  
  58. [target]
  59. ##use a tab-delimited .BED file to specify capture regions (control dataset is needed to use this option):
  60. #captureRegions = /bioinfo/users/vboeva/Desktop/testChr19/capture.bed

后面的target标签应该是做靶向捕获数据或者只想分析某个区域数据时才需要的，而BAF标签则是需要分析LOH等内容时才需要。要进行BAF分析时，输入文件格式必须时pileup而不是bam。

general标签中的chrLenFile需要的是染色体的长度，可以直接使用samtools faidx 生成的fai索引（未测试），官网也提供了hg19的这个文件。如果是使用GATK的b37参考基因组，则把chr去掉即可。

chrFiles需要填的是fasta格式文件的路径，建议是该路径下仅包含reference的各个染色体的fasta文件，不要包含其他基因组的fasta文件，命名要与chrLenFile中的染色体名相同。如chr1.fa、chr2.fa等。

可以使用samtools对参考基因组进行提取，把每个染色体拆分出来

  1. samtools faidx hg19.fa chr1 > seq_chrom/chr1.fa

GCcontentProfile则是需要划分的每个bin的GC比例用于校正，获得的方法是使用gccount软件，conf.txt文件与freec的用同一个即可。

  1. gccount -conf conf.txt

后续分析

软件的速度还是比较快的，上述分析完成后，会生成后缀为_CNVs、_ratio.txt、_sample.cpn、_info.txt等文件。结果在CNVs文件中。

再使用assess_significance.R脚本补充p值，这个脚本需要安装R包“rtracklayer”

  1. cat FREEC-11.6/scripts/assess_significance.R | R --slave --args test.bam_CNVs test.bam_ratio.txt

可使用makeGraph.R画图

  1. cat FREEC-11.6/scripts/makeGraph.R | R --slave --args 2 test.bam_ratio.txt

Control-FREEC的使用

将ratio文件转为bed格式或circos格式

  1. perl FREEC-11.6/scripts/freec2bed.pl -f test.bam_ratio.txt > test.bam_ratio.bed
  2. perl FREEC-11.6/scripts/freec2circos.pl -f test.bam_ratio.txt > test.bam_ratio.circos

发表评论

匿名网友

拖动滑块以完成验证
加载失败