runAbsoluteCN {PureCN} | R Documentation

## Run PureCN implementation of ABSOLUTE

### Description

This function takes as input tumor and normal control coverage data and a VCF containing allelic fractions of germline variants and somatic mutations. Normal control does not need to be from the same patient. In case VCF does not contain somatic status, it should contain dbSNP and optionally COSMIC annotation. Returns purity and ploidy combinations, sorted by likelihood score. Provides copy number and LOH data, by both gene and genomic region.

### Usage

    runAbsoluteCN(normal.coverage.file = NULL, tumor.coverage.file = NULL,
      log.ratio = NULL, seg.file = NULL, seg.file.sdev = 0.4,
      vcf.file = NULL, normalDB = NULL, genome, centromeres = NULL,
      sex = c("?", "F", "M", "diploid"), fun.filterVcf = filterVcfMuTect,
      args.filterVcf = list(), fun.setPriorVcf = setPriorVcf,
      args.setPriorVcf = list(), fun.setMappingBiasVcf = setMappingBiasVcf,
      args.setMappingBiasVcf = list(), fun.filterTargets = filterTargets,
      args.filterTargets = list(), fun.segmentation = segmentationCBS,
      args.segmentation = list(), fun.focal = findFocal, args.focal = list(),
      sampleid = NULL, min.ploidy = 1, max.ploidy = 6, test.num.copy = 0:7,
      test.purity = seq(0.15, 0.95, by = 0.01), prior.purity = NULL,
      prior.K = 0.999, prior.contamination = 0.01,
      max.candidate.solutions = 20, candidates = NULL, min.coverage = 15,
      max.coverage.vcf = 300, max.non.clonal = 0.2,
      max.homozygous.loss = c(0.05, 1e+07), non.clonal.M = 1/3,
      max.mapping.bias = 0.8, max.pon = 3, iterations = 30,
      min.variants.segment = 5, log.ratio.calibration = 0.1,
      smooth.log.ratio = TRUE, remove.off.target.snvs = NULL,
      model.homozygous = FALSE, error = 0.001, gc.gene.file = NULL,
      max.dropout = c(0.95, 1.1), min.logr.sdev = 0.15, max.logr.sdev = 0.75,
      max.segments = 300, min.gof = 0.8, plot.cnv = TRUE,
      cosmic.vcf.file = NULL, model = c("beta", "betabin"),
      post.optimize = FALSE, log.file = NULL, verbose = TRUE)

### Arguments

`normal.coverage.file` | Coverage file of normal control (optional if log.ratio is provided - then it will be only used to filter low coverage exons). Should be already GC-normalized with `correctCoverageBias`. Needs to be either a file name or data read with the `readCoverageFile` function.
`tumor.coverage.file` | Coverage file of tumor. If `NULL`, requires `seg.file` and an interval file via `gc.gene.file`. Should be already GC-normalized with `correctCoverageBias`. Needs to be either a file name or data read with the `readCoverageFile` function.
`log.ratio` | Copy number log-ratios for all exons in the coverage files. If `NULL`, calculated based on coverage files.
`seg.file` | Segmented data. Optional, to support matched SNP6 data. If `NULL`, use coverage files or `log.ratio` to segment the data.
`seg.file.sdev` | If `seg.file` provided, the log-ratio standard deviation, used to model likelihood of sub-clonal copy number events.
`vcf.file` | VCF file, tested with ‘MuTect 1’ output files. Optional, but typically needed to select between local optima of similar likelihood. Can also be a `CollapsedVCF`, read with the `readVcf` function. Requires a DB info flag for dbSNP membership. The default `fun.setPriorVcf` function will also look for a Cosmic.CNT slot (see `cosmic.vcf.file`), containing the hits in the COSMIC database. Again, do not expect very useful results without a VCF file.
`normalDB` | Normal database, created with `createNormalDatabase`. If provided, used to calculate gene-level p-values (requires `Gene` column in `gc.gene.file`) and to filter targets with low coverage in the pool of normal samples.
`genome` | Genome version, for example hg19.
`centromeres` | A `GRanges` object with centromere positions. If `NULL`, use pre-stored positions for genome versions hg18, hg19 and hg38.
`sex` | Sex of sample. If `?`, detect using `getSexFromCoverage` function and default parameters. Default parameters might not work well with every assay and might need to be tuned. If set to diploid, then PureCN will assume all chromosomes are diploid and will not try to detect sex.
`fun.filterVcf` | Function for filtering variants. Expected output is a list with elements `vcf` (`CollapsedVCF`), flag (`logical(1)`) and `flag_comment` (`character(1)`). The flags will be added to the output data and can be used to warn users, for example when samples look too noisy. Default filter will remove variants flagged by MuTect, but will keep germline variants. If ran in matched normal mode, it will by default use somatic status of variants and filter non-somatic calls with allelic fraction significantly different from 0.5 in normal. Defaults to `filterVcfMuTect`, which in turn also calls `filterVcfBasic`.
`args.filterVcf` | Arguments for variant filtering function. Arguments `vcf`, `tumor.id.in.vcf`, `min.coverage`, `model.homozygous` and `error` are required in the filter function and are automatically set.
`fun.setPriorVcf` | Function to set prior for somatic status for each variant in the VCF. Defaults to `setPriorVcf`.
`args.setPriorVcf` | Arguments for somatic prior function.
`fun.setMappingBiasVcf` | Function to set mapping bias for each variant in the VCF. Defaults to `setMappingBiasVcf`.
`args.setMappingBiasVcf` | Arguments for mapping bias function.
`fun.filterTargets` | Function for filtering low-quality targets in the coverage files. Needs to return a `logical` vector whether an interval should be used for segmentation. Defaults to `filterTargets`.
`args.filterTargets` | Arguments for target filtering function. Arguments `normal`, `tumor`, `log.ratio`, `min.coverage` `seg.file` and `normalDB` are required and automatically set.
`fun.segmentation` | Function for segmenting the copy number log-ratios. Expected return value is a `data.frame` representation of the segmentation. Defaults to `segmentationCBS`.
`args.segmentation` | Arguments for segmentation function. Arguments `normal`, `tumor`, `log.ratio`, `plot.cnv`, `sampleid`, `vcf`, `tumor.id.in.vcf`, `centromeres` are required in the segmentation function and automatically set.
`fun.focal` | Function for identifying focal amplifications. Defaults to `findFocal`.
`args.focal` | Arguments for focal amplification function.
`sampleid` | Sample id, provided in output files etc.
`min.ploidy` | Minimum ploidy to be considered.
`max.ploidy` | Maximum ploidy to be considered.
`test.num.copy` | Copy numbers tested in the grid search. Note that focal amplifications can have much higher copy numbers, but they will be labeled as subclonal (because they do not fit the integer copy numbers).
`test.purity` | Considered tumor purity values.
`prior.purity` | `numeric(length(test.purity))` with priors for tested purity values. If `NULL`, use flat priors.
`prior.K` | This defines the prior probability that the multiplicity of a SNV corresponds to either the maternal or the paternal copy number (for somatic variants additionally to a multiplicity of 1). For perfect segmentations, this value would be 1; values smaller than 1 thus may provide some robustness against segmentation errors.
`prior.contamination` | The prior probability that a known SNP is from a different individual.
`max.candidate.solutions` | Number of local optima considered in optimization and variant fitting steps. If there are too many local optima, it will use specified number of top candidate solutions, but will also include all optima close to diploid, because silent genomes have often lots of local optima.
`candidates` | Candidates to optimize from a previous run (`return.object$candidates`). If `NULL`, do 2D grid search and find local optima.
`min.coverage` | Minimum coverage in both normal and tumor. Targets and variants with lower coverage are ignored. This value is provided to the `args.filterTargets` and `args.filterVcf` lists, but can be overwritten in these lists if different cutoffs for the coverage and variant filters are wanted. To increase the sensitivity of homozygous deletions in high purity samples, the coverage cutoff in tumor is automatically lowered by 50 percent if the normal coverage is high.
`max.coverage.vcf` | This will set the maximum number of reads in the SNV fitting. This is to avoid that small non-reference biases that come apparent only at high coverages have a dramatic influence on likelihood scores.
`max.non.clonal` | Maximum genomic fraction assigned to a subclonal copy number state.
`max.homozygous.loss` | `double(2)` with maximum genomic fraction assigned to homozygous loss and maximum size of a homozygous loss segment. These are set to a fairly high default value to not exclude correct solutions, especially in noisy segmentations.
`non.clonal.M` | Average expected cellular fraction of sub-clonal somatic mutations. This is to calculate expected allelic fractions of a single sub-clonal bin for variants. For all somatic variants, more accurate cellular fractions are calculated.
`max.mapping.bias` | Exclude variants with high mapping bias from the likelihood score calculation. Note that bias is reported on an inverse scale; a variant with mapping bias of 1 has no bias.
`max.pon` | Exclude variants found more than `max.pon` times in pool of normals and not in dbSNP. Requires `normal.panel.vcf.file` in `setMappingBiasVcf`. Should be set to a value high enough to be much more likely an artifact and not a true germline variant not present in dbSNP.
`iterations` | Maximum number of iterations in the Simulated Annealing copy number fit optimization. Note that this an integer optimization problem that should converge quickly. Allowed range is 10 to 250.
`min.variants.segment` | Flag segments with fewer variants. The minor copy number estimation is not reliable with insufficient variants.
`log.ratio.calibration` | Re-calibrate log-ratios in the window `sd(log.ratio)*log.ratio.calibration`.
`smooth.log.ratio` | Smooth `log.ratio` using the `DNAcopy` package.
`remove.off.target.snvs` | Deprecated. Use the corresponding argument in `args.filterVcf`.
`model.homozygous` | Homozygous germline SNPs are uninformative and by default removed. In 100 percent pure samples such as cell lines, however, heterozygous germline SNPs appear homozygous in case of LOH. Setting this parameter to `TRUE` will keep homozygous SNPs and include a homozygous SNP state in the likelihood model. Not necessary when matched normal samples are available.
`error` | Estimated sequencing error rate. Used to calculate minimum number of supporting reads for variants using `calculatePowerDetectSomatic`. Also used to calculate the probability of homozygous SNP allelic fractions (assuming reference reads are sequencing errors).
`gc.gene.file` | A mapping file that assigns GC content and gene symbols to each exon in the coverage files. Used for generating gene-level calls. First column in format CHR:START-END. Second column GC content (0 to 1). Third column gene symbol. This file can be generated with the ‘GATK GCContentByInterval’ tool or with the `calculateGCContentByInterval` function.
`max.dropout` | Measures GC bias as ratio of coverage in AT-rich (GC < 0.5) versus GC-rich on-target regions (GC >= 0.5). High drop-out might indicate that data was not GC-normalized or that the sample quality might be insufficient. Requires `gc.gene.file`.
`min.logr.sdev` | Minimum log-ratio standard deviation used in the model. Useful to make fitting more robust to outliers in very clean data.
`max.logr.sdev` | Flag noisy samples with segment log-ratio standard deviation larger than this. Assay specific and needs to be calibrated.
`max.segments` | Flag noisy samples with a large number of segments. Assay specific and needs to be calibrated.
`min.gof` | Flag purity/ploidy solutions with poor fit.
`plot.cnv` | Generate segmentation plots.
`cosmic.vcf.file` | Add a `Cosmic.CNT` info field to the provided `vcf.file` using a VCF file containing the COSMIC database. The default `fun.setPriorVcf` function will give variants found in the COSMIC database a higher prior probability of being somatic. Not used in likelhood model when matched normal is available in `vcf.file`. Should be compressed and indexed with bgzip and tabix, respectively.
`model` | Use either a beta or a beta-binomial distribution for fitting observed to expected allelic fractions of alterations in `vcf.file`. The latter can be useful to account for significant overdispersion, for example due to mapping biases when no pool of normals is available or due to other unmodeled biases, e.g. amplification biases. The amount of expected overdispersion can be controlled via the `max.coverage.vcf` argument (the higher, the less expected bias).
`post.optimize` | Optimize purity using final SCNA-fit and variants. This might take a long time when lots of variants need to be fitted, but will typically result in a slightly more accurate purity, especially for rather silent genomes or very low purities. Otherwise, it will just use the purity determined via the SCNA-fit.
`log.file` | If not `NULL`, store verbose output to file.
`verbose` | Verbose output.

### Value

A list with elements

`candidates` | Results of the grid search.
`results` | All local optima, sorted by final rank.
`input` | The input data.

### Author(s)

Markus Riester

### References

Riester et al. (2016). PureCN: Copy number calling and SNV classification using targeted short read sequencing. Source Code for Biology and Medicine, 11, pp. 13.

Carter et al. (2012), Absolute quantification of somatic DNA alterations in human cancer. Nature Biotechnology.

### See Also

`correctCoverageBias` `segmentationCBS` `calculatePowerDetectSomatic`

### Examples

    normal.coverage.file <- system.file('extdata', 'example_normal.txt',
        package='PureCN')
    tumor.coverage.file <- system.file('extdata', 'example_tumor.txt',
        package='PureCN')
    vcf.file <- system.file('extdata', 'example_vcf.vcf',
        package='PureCN')
    gc.gene.file <- system.file('extdata', 'example_gc.gene.file.txt',
        package='PureCN')

    # The max.candidate.solutions, max.ploidy and test.purity parameters are set to
    # non-default values to speed-up this example.  This is not a good idea for real
    # samples.
    ret <-runAbsoluteCN(normal.coverage.file=normal.coverage.file,
        tumor.coverage.file=tumor.coverage.file, genome='hg19', vcf.file=vcf.file,
        sampleid='Sample1', gc.gene.file=gc.gene.file,
        max.ploidy=4, test.purity=seq(0.3,0.7,by=0.05), max.candidate.solutions=1)

    # If a high-quality segmentation was obtained with third-party tools:
    seg.file <- system.file('extdata', 'example_seg.txt',
        package = 'PureCN')

    # By default, PureCN will re-segment the data, for example to identify
    # regions of copy number neutral LOH. If this is not wanted, we can provide
    # a minimal segmentation function which just returns the provided one:
    funSeg <- function(seg, ...) return(seg)

    res <- runAbsoluteCN(seg.file=seg.file, fun.segmentation=funSeg, max.ploidy = 4,
        test.purity = seq(0.3, 0.7, by = 0.05), max.candidate.solutions=1,
        genome='hg19', gc.gene.file=gc.gene.file)

[Package PureCN version 1.8.1 Index]