Skip to content

ShixiangWang/gcap

Repository files navigation

GCAP: Gene-level Circular Amplicon Prediction

gcap status badge lifecycle Badge

In a nutshell, gcap provides an end-to-end workflow for predicting circular amplicon (also known as ecDNA, extra-chromosomal DNA ) in gene level with machine learning approach, then classifying cancer samples into different focal amplification (fCNA) types, based on input from WES (tumor-normal paired BAM, with corresponding .bai index files) data, allele specific copy number data (e.g., results from ASCAT or Sequenza), or even absolute integer copy number data (e.g., results from ABSOLUTE). The former two data sources are preferred as input of gcap .

Installation

Install alleleCount (WES bam data only)

alleleCount is required to run ASCAT on WES bam data, if you haven't installed conda or miniconda, please install firstly, then install the alleleCount in terminal with:

conda create -n cancerit -c bioconda cancerit-allelecount

NOTE: gcap set the default alleleCount as the ~/miniconda3/envs/cancerit/bin/alleleCounter, if you use conda or other approaches, please set the path when you use corresponding functions.

Install ASCAT (required)

Latest ASCAT v3

From v1.2, GCAP uses the latest version of ASCAT. Install ASCAT v3 in R console from GitHub with:

# install.packages("remotes")
remotes::install_github('VanLoo-lab/ascat/ASCAT')

We have provided generated reference files above, but sometimes you may want to generate the reference data for yourself, in such case, please refer to https://github.com/VanLoo-lab/ascat for generating the required allele-specific copy number data.

Reference files:

The reference files are required in ASCAT for copy number calling.

The prediction model was built with data on the top of hg38 genome build, so hg38-based BAM file input is more recommended.

A fixed version of ASCAT v3

In our manuscript, we used a fixed version of ASCAT for the GCAP data pre-processing (modified and adapted for GCAP workflow in HPC). It does not fit the R version >=4.3.

# This is a forked version ASCAT
remotes::install_github("ShixiangWang/ascat@v3-for-gcap-v1", subdir = "ASCAT")
# A ASCAT version with loose SAM flag, useful sometimes
# remotes::install_github("ShixiangWang/ascat@v3-f1", subdir = "ASCAT")
# See https://github.com/ShixiangWang/gcap/issues/27

Reference files:

Alternatives to ASCAT

For the latest version of GCAP, sequenza or facets are supported for preprocessing the bam data, please refer to the provided links for usage.

Install GCAP (required)

Install gcap in R console:

# r-universe
install.packages('gcap', repos = c('https://shixiangwang.r-universe.dev', 'https://cloud.r-project.org'))

# or GitHub
remotes::install_github("ShixiangWang/gcap")

To work with the fixed version of ASCAT, you have to install version commits no more newer than 42f216d (tag v1.1.5), i.e., please use remotes::install_github("ShixiangWang/gcap@v1.1.5") and the R version should below v4.3.

If you would like to use CLI program in Shell terminal, run the following code in your R console after installation:

gcap::deploy()

Two scripts gcap-bam.R and gcap-ascn.R shall be linked to your path /usr/local/bin/. You can use one of them based on you input data.

NOTE

For users with package GetoptLong version >= 1.1.0, a main command is implemented and also linked to /usr/local/bin/ when calling deploy(). So you can type gcap as a unified interface.

$ gcap
gcap (v1.0.0)
Usage: gcap [command] [options]

Commands:
  bam     Run GCAP workflow with tumor-normal paired BAM files
  ascn    Run GCAP workflow with curated allele-specific copy number data

----------
Citation:
  GCAP
URL:
  https://github.com/ShixiangWang/gcap

NOTE: gcap use XGBOOST < 1.6, if you have installed a latest version, you can install the specified version with:

install.packages("https://cran.r-project.org/src/contrib/Archive/xgboost/xgboost_1.5.2.1.tar.gz", repos = NULL)

Example

Run the following code to see a quick example:

library(gcap)

data("ascn")
rv <- gcap.ASCNworkflow(ascn, outdir = tempdir(), model = "XGB11")
rv

Pipeline (WES bam data only)

To run gcap from bam files, a machine with at least 80GB RAM is required for the allelecount process. If you set multiple threads, please note the parallel computation is used in part of the workflow. You should balance the nthread setting and the computing power your machine provides by yourself.

It generally takes ~0.5h to finish one case (tumor-normal pair).

In our practice, when we want to process multiple cases, set nthread = 22 and directly let gcap handle multiple cases (instead of writing a loop yourself) is good enough.

A recommended setting for Slurm is given as:

#!/bin/bash
#SBATCH -N 1
#SBATCH -o output-%J.o
#SBATCH -n 22
#SBATCH --mem=102400

Templates of practical calling command with provided hg38 and hg19 annotations are given below:

# hg38 ----------------
gcap.workflow(
  tumourseqfile = tfile, normalseqfile = nfile, tumourname = tn, normalname = nn, jobname = id,
  outdir = outdir,
  allelecounter_exe = "~/miniconda3/envs/cancerit/bin/alleleCounter", 
  g1000allelesprefix = file.path(
    "/data/wsx/data/1000G_loci_hg38/",
    "1kg.phase3.v5a_GRCh38nounref_allele_index_chr"
  ), 
  g1000lociprefix = file.path("/data/wsx/data/1000G_loci_hg38/",
                              "1kg.phase3.v5a_GRCh38nounref_loci_chrstring_chr"
  ),
  GCcontentfile = "/data/wsx/data/GC_correction_hg38.txt",
  replictimingfile = "/data/wsx/data/RT_correction_hg38.txt",
  skip_finished_ASCAT = TRUE,
  skip_ascat_call = FALSE,
  result_file_prefix = "xxx",
  extra_info = df,
  include_type = FALSE,
  genome_build = "hg38",
  model = "XGB11"
)

# hg19 ----------------
gcap.workflow(
  tumourseqfile = tfile, normalseqfile = nfile, tumourname = tn, normalname = nn, jobname = id,
  outdir = outdir,
  allelecounter_exe = "~/miniconda3/envs/cancerit/bin/alleleCounter", g1000allelesprefix = file.path(
    "/data/wsx/data/1000G_loci_hg19/",
    "1000genomesAlleles2012_chr"
  ), g1000lociprefix = file.path("/data/wsx/data/1000G_loci_hg19/", "1000genomesloci2012chrstring_chr"),
  GCcontentfile = "/data/wsx/data/GC_correction_hg19.txt", replictimingfile = "/data/wsx/data/RT_correction_hg19.txt",
  skip_finished_ASCAT = TRUE,
  skip_ascat_call = FALSE,
  result_file_prefix = "xxx",
  extra_info = NULL,
  include_type = FALSE,
  genome_build = "hg19",
  model = "XGB11"
)

Pipeline (Allele specific or absolute copy number data)

Please refer to ?gcap.ASCNworkflow().

Functions

For more custom and advanced control of the analysis, you can read the structured documentation at package site.

Logging

For better debugging and rechecking. The logging information of your operation with gcap would be saved into an independent file. You can use the following commands to get the file path and print logging message. Please note you have to use ::: to access these functions as they are not exported from gcap.

> gcap:::get_log_file()
[1] "~/Library/Logs/gcap/gcap.log"
> gcap:::cat_log_file()

Docker image

A docker image is available in ghcr along with its corresponding Dockerfile. This image comes pre-installed with all the necessary software (including R, alleleCounter, ASCAT, sequenza-utils, samtools, etc.) with conda environment pre-configured — just docker run and the gcap command works out of the box.

However, users are responsible for mapping the required reference files and input data files on their own. The Dockerfile can be customized according to the user's specific requirements, as permitted by the license we provide.

docker pull ghcr.io/shixiangwang/gcap:latest

# Quick test
docker run ghcr.io/shixiangwang/gcap:latest --help

# Run with mounted data
docker run -v /path/to/data:/data ghcr.io/shixiangwang/gcap:latest bam [options]

Note: In the Docker container, alleleCounter is installed at /opt/conda/bin/alleleCounter. Set allelecounter_exe = "/opt/conda/bin/alleleCounter" when using gcap.workflow() inside the container.

Related tools

  • gcaputils: gcap utils for downstream analysis and visualization.
  • DoAbsolute: Automate Absolute Copy Number Calling using 'ABSOLUTE' package.
  • sigminer: An easy-to-use and scalable toolkit for genomic alteration signature (a.k.a. mutational signature) analysis and visualization in R.

Output

gcap workflows (gcap.workflow(), gcap.ASCNworkflow()) return an fCNA R6 object containing prediction results at both gene and sample level. Two CSV files are also saved to disk.

Understanding the Results

The fCNA object provides two views of your data:

Field Description
$data Gene-level predictions: one row per gene per sample. Key columns: sample, gene_id, band (cytoband), total_cn, minor_cn, prob (ecDNA probability), gene_class ("circular" or "noncircular")
$sample_summary Sample-level classification: one row per sample. Key columns: sample, purity, ploidy, class ("circular", "noncircular", or "nofocal"), ec_genes (count of circular genes), ec_cytobands

Sample classification logic: A sample is classified as "circular" if at least one cytoband has a high-confidence circular amplicon (default probability threshold ≥ 0.6). Multiple genes within the same cytoband are aggregated — the sample is circular if any cytoband passes the threshold. "nofocal" means no focal amplification was detected at all.

Exploring Results

rv <- gcap.ASCNworkflow(ascn, outdir = tempdir(), model = "XGB11")

# Print summary
rv                    # shows sample & gene counts, sample classification breakdown

# Get gene-level data
head(rv$data)         # all predicted genes with probabilities and classifications
table(rv$data$gene_class)  # count of circular vs noncircular genes

# Get sample-level summary
head(rv$sample_summary)    # one row per sample
table(rv$sample_summary$class)  # count of circular/noncircular/nofocal samples

# Which samples have ecDNA (circular)?
rv$sample_summary[sample_summary$class == "circular", ]

# What genes are predicted circular across all samples?
rv$getGeneSummary()         # frequency table of gene classifications
rv$getCytobandSummary()     # frequency table of cytoband classifications

# Filter to specific samples
subset_rv <- rv$subset(sample %in% c("TCGA-02-2485-01", "TCGA-05-4395-01"))
subset_rv  # re-prints summary for subset

# Filter to high-confidence circular genes only
high_conf <- rv$subset(total_cn > 10)  # genes with high-level amplification

Saved Files

File Content
{prefix}_prediction_result.rds Raw prediction data (RDS format)
{prefix}_fCNA_records.csv Gene-level predictions table
{prefix}_sample_info.csv Sample-level classification summary

Citations

Wang, S., Wu, C. Y., He, M. M., Yong, J. X., Chen, Y. X., Qian, L. M., ... & Zhao, Q. (2024). Machine learning-based extrachromosomal DNA identification in large-scale cohorts reveals its clinical implications in cancer. Nature Communications, 15(1), 1-17. https://doi.org/10.1038/s41467-024-45479-6

LICENSE

This software and associated documentation files (the “Software”) are protected by copyright. This Software is provided “as is” (at your own risk) for internal non-commercial academic research purposes only. Please read the Non-Commercial Academic License in detail before downloading a copy. By installing or using this Software, you agree to be bound by the terms and conditions of the Non-Commercial Academic License.

All commercial use of the Software or any modification, manipulation or derivative of the Software, including but not limited to transfer, sale or licence to a commercial third party or use on behalf of a commercial third party (including but not limited to use as part of a service supplied to any third party for financial reward) is strictly prohibited and requires a commercial use licence.

Patent Protection

This software is protected by the P. R. China patent 202211067952.6. The following intellectual property is protected under patent law and is NOT covered by the academic non-commercial license — a separate commercial license is required for any use beyond personal academic research:

  • Modeling methods: the machine learning pipeline for gene-level circular amplicon (ecDNA) prediction, feature engineering strategies, and focal amplification (fCNA) classification algorithms.
  • Modeling features and approaches: derived genomic features, copy number variation metrics, and their computational extraction methods.
  • Model construction: XGBoost-based model architecture, training methodology, hyperparameter configurations, and ensemble strategies.
  • Model files: pre-trained model objects, serialized model binaries, model parameters, and any derivative works thereof.

For further information or to obtain a commercial license, please email wangshx@csu.edu.cn or zhaoqi@sysucc.org.cn.

About

GCAP (Gene-level Circular Amplicon Prediction) firstly implements extrachromosomal DNA detection from whole-exome-sequencing (WES) data and absolute copy number profiles. https://shixiangwang.r-universe.dev/gcap

Topics

Resources

License

Stars

21 stars

Watchers

2 watching

Forks

Packages

 
 
 

Contributors