The example or result data in the result folder can be downloaded here. The example or result data in the result_merge folder can be downloaded here.

Calculate the number of fragments with copy number amplification in each cell

We do not provide example data for the estimate_Copy_num function because the example data for a single sample exceeds 50GB. The input file structure required to run this function for a single sample is as follows:

## eccDNAscope_example/sample_data/AK1
## └── result_cellranger-atac
##     ├── AK1.log
##     ├── AK1_fragments.tsv.gz
##     ├── AK1_fragments.tsv.gz.tbi
##     ├── analysis
##     ├── cloupe.cloupe
##     ├── cut_sites.bigwig
##     ├── filtered_peak_bc_matrix
##     ├── filtered_peak_bc_matrix.h5
##     ├── filtered_tf_bc_matrix
##     ├── filtered_tf_bc_matrix.h5
##     ├── peak_annotation.tsv
##     ├── peak_motif_mapping.bed
##     ├── peaks.bed
##     ├── possorted_bam.bam
##     ├── possorted_bam.bam.bai
##     ├── possorted_bam.md5
##     ├── raw_peak_bc_matrix
##     ├── raw_peak_bc_matrix.h5
##     ├── singlecell.csv
##     ├── summary.csv
##     ├── summary.json
##     └── web_summary.html

Attention: eccDNAscope_example/sample_data/AK1 is the storage path for a single sample, while result_cellranger-atac is the relative path under AK1 where the Cell Ranger results are stored.The cellranger directory must be generated using the standard Cell Ranger-ATAC pipeline. The following files are required to run this function: singlecell.csv, possorted_bam.bam, possorted_bam.bam.bai

library(parallel)
samples <- c(
  "eccDNAscope_example/sample_data/AK1",
  "eccDNAscope_example/sample_data/CSCC1")

# Parallelly calculate the CNV_copy feature for all samples using the estimate_Copy_num function.
results <- mclapply(samples, function(dir) {
  estimate_Copy_num(
    eccDNA.min.copy = 5,
    sample.dir = dir,
    cellranger.result.dir="result_cellranger-atac",
    blacklist.file="eccDNAscope_example/reference/blacklist_repeats_segdups_rmsk_hg38.bed",
    bedtools.dir=".conda/envs/eccDNAscope/bin/bedtools",
    output.dir="eccDNAscope_example/result"
  )
}, mc.cores = length(samples))

# Merge results from all samples.
CNV_copy <- do.call(c, lapply(results, function(x) x$CNV_copy))
CNV_allcopy <- do.call(c, lapply(results, function(x) x$CNV_allcopy))
save(CNV_copy,CNV_allcopy,file='eccDNAscope_example/result_merge/CNV_statistic.Rdata')

Assign CNV_copy information to the eccDNA object

load("eccDNAscope_example/result_merge/CNV_statistic.Rdata")
combined_atac <- readRDS("eccDNAscope_example/result_merge/combined_atac_annote.rds")
combined_atac <- assign_eccDNA_features(eccDNA.object=combined_atac,
                                       cell_type_col="celltype_ATAC",
                                       feature_name = "CNV_copy")

Visualize the distribution of eccDNA features across cells.

# CNV_copy distribution in all cells
combined_atac@misc$p_CNV_copy_distribution

# CNV_copy distribution in each cell type
combined_atac@misc$p_CNV_copy_celltype


Calculate the number of adapter(breakpoint) of eccDNA in each cell

Before running this function, the estimate_Copy_num function must have been executed with the same output.dir parameter.

library(parallel)
samples <- c("eccDNAscope_example/sample_data/AK1",
             "eccDNAscope_example/sample_data/CSCC1")

# Parallelly calculate the adapter_num feature for all samples using the estimate_Adapter_num function.
results <- mclapply(samples, function(dir) {
  estimate_Adapter_num(
    sample.dir = dir,
    cellranger.result.dir="result_cellranger-atac",
    samtools.dir=".conda/envs/eccDNAscope/bin/samtools",
    bedtools.dir=".conda/envs/eccDNAscope/bin/bedtools",
    output.dir="eccDNAscope_example/result"  # The value of the output.dir parameter must be the same as the output.dir parameter used in the estimate_Copy_num function.
  )
}, mc.cores = length(samples))

# Merge results from all samples.
adapter_num <- do.call(c, results)
save(adapter_num,file='eccDNAscope_example/result_merge/adapter_num.Rdata')

Assign adapter_num information to the eccDNA object

load("eccDNAscope_example/result_merge/adapter_num.Rdata")
combined_atac <- assign_eccDNA_features(eccDNA.object=combined_atac,
                                       cell_type_col="celltype_ATAC",
                                       feature_name = "adapter_num")

Visualize the distribution of eccDNA features across cells.

# adapter_num distribution in all cells
combined_atac@misc$p_adapter_num_distribution

# adapter_num distribution in each cell type
combined_atac@misc$p_adapter_num_celltype


Calculate the number of eccDNA genes in each cell

library(parallel)
sample.file.names <- c("AK1","CSCC1")
# Parallelly calculate the eccDNA_gene_num feature for all samples using the estimate_eccDNA_gene_num function.
results <- mclapply(sample.file.names, function(sample) {
estimate_eccDNA_gene_num(
   sample.file.name = sample,
   adapter.result.dir="eccDNAscope_example/result"
 )
}, mc.cores = length(sample.file.names))

# Merge results from all samples.
eccDNA_gene_num <- do.call(c, results)
save(eccDNA_gene_num,file='eccDNAscope_example/result_merge/eccDNA_gene_num.Rdata')

Assign eccDNA_gene_num information to the eccDNA object

load("eccDNAscope_example/result_merge/eccDNA_gene_num.Rdata")
combined_atac <- assign_eccDNA_features(eccDNA.object=combined_atac,
                                       cell_type_col="celltype_ATAC",
                                       feature_name = "eccDNA_gene_num")

Visualize the distribution of eccDNA features across cells.

# eccDNA_gene_num distribution in all cells
combined_atac@misc$p_eccDNA_gene_num_distribution

# eccDNA_gene_num distribution in each cell type
combined_atac@misc$p_eccDNA_gene_num_celltype


Plot eccDNA structure

Select cells eligible for eccDNA visualization from a specific sample.

cell_type <- combined_atac$celltype_ATAC
names(cell_type) <- rownames(combined_atac@meta.data)
load("eccDNAscope_example/result/CSCC1/final_copy.Rdata") # If you have executed the first three steps, the file will be automatically generated in the output path.
final_copy_filter <- get_eccDNA_cells(final_copy=final_copy,cell_type=cell_type)

# The cell names you can choose to plot
new_exist_cell<-c()
for(i in 1:dim(final_copy_filter)[1]){
  temp_data<-final_copy_filter[i,10:13]
  if(length(grep('-1$',temp_data))!=0){
    temp_index<-strsplit(temp_data[grep('-1$',temp_data)],':')
    new_exist_cell[i]<-temp_index[[1]][[3]]
  }else{
    new_exist_cell[i]<-0
  }
}
sort(table(new_exist_cell),decreasing = T)[1:5]
## TAGCACAAGGTTCTCA-1               <NA>               <NA>               <NA> 
##                  4                 NA                 NA                 NA 
##               <NA> 
##                 NA

Visualizing eccDNA in specific cell

p_one_cell <- plot_eccDNAs_in_one_cell(cell.name="TAGCACAAGGTTCTCA-1",final_copy_filter=final_copy_filter)
## >> preparing features information...      2025-8-25 21:49:42 
## >> identifying nearest features...        2025-8-25 21:49:43 
## >> calculating distance from peak to TSS...   2025-8-25 21:49:43 
## >> assigning genomic annotation...        2025-8-25 21:49:43 
## >> adding gene annotation...          2025-8-25 21:50:07 
## >> assigning chromosome lengths           2025-8-25 21:50:07 
## >> done...                    2025-8-25 21:50:07

Filter out eccDNA-related genes eligible in the specific cell.

table(p_one_cell$one_cell$SYMBOL)
## 
## PTMA 
##    4

Plot the amplification structure of a single eccDNA-related gene in a specific cell.

p <- plot_single_eccDNA_in_one_cell(gene.name="PTMA",one_cell=p_one_cell$one_cell)


Visualize the genomic locations of eccDNA

Integrate eccDNA-related gene amplification information across all samples.

cell_type <- combined_atac$celltype_ATAC
names(cell_type) <- rownames(combined_atac@meta.data)
all_annotation <- integrate_eccDNA_amplified_genes(adapter.result.dir="eccDNAscope_example/result",
                                                                     cell_type=cell_type)

Plot the relationship between breakpoints and copy numbers of eccDNA-amplified genes

plot_eccDNA_amplified_genes(all_annotation=all_annotation)

Plot the distribution of eccDNA genomic regions

p <- plot_eccDNA_genomic_distribution(all_annotation=all_annotation)
## >> preparing features information...      2025-8-25 21:50:12 
## >> identifying nearest features...        2025-8-25 21:50:12 
## >> calculating distance from peak to TSS...   2025-8-25 21:50:13 
## >> assigning genomic annotation...        2025-8-25 21:50:13 
## >> adding gene annotation...          2025-8-25 21:50:18
## 'select()' returned 1:many mapping between keys and columns
## >> assigning chromosome lengths           2025-8-25 21:50:18 
## >> done...                    2025-8-25 21:50:18

p$annotation_res