Using ExomeDepth for GRCH38 processed samples to call CNVs

Question

[nkausthu]

Hello,

I use ExomeDepth for CNV analysis for GRCh37 processed samples and now the processing pipeline is shifted to GRCh38. Unfortunately I couldn't find any protocol to use ExomeDepth for GRCh38. I would really appreciate if someone can help me with this. Looking forward for the solutions. Thank you.

[Jorge Amigo]

Although someone may have decided to solve this issue and hopefully would like to share his experience here, I just wanted to let you know that the author unfortunately hasn't answered such question for almost a year now.

Answer 1

[Joakim]

The only difference would be the annotations, instead of using bedframes from data(genes.hg19) and data(exons.hg19) in ExomeDepth, I got them from the UCSC Table Browser for hg38 (http://genome.ucsc.edu/cgi-bin/hgTables).

The only info they contain are: chromosome start end name

..and then run as before. Change bed.frame = exons.hg19 to the exon data you downloaded.

I called the exon data "geneid":

my.counts <- getBamCounts(bed.frame = geneid,
                          bam.files = my.bam,
                          include.chr = F)

In the end you can annotate also with genes from UCSC:

Annotate Genes

genes.hg38 <- read.csv("Genes.csv", header=TRUE, sep = ";")

genes.hg38.GRanges <- GenomicRanges::GRanges(seqnames = genes.hg38$chromosome,
IRanges::IRanges(start=genes.hg38$start,end=genes.hg38$end),
                                             names = genes.hg38$geneName)

all.exons <- AnnotateExtra(x = all.exons,
                           reference.annotation = genes.hg38.GRanges,
                           min.overlap = 0.1,
                           column.name = 'genes.hg38')

[Joel Wallenius]

Hey, I'm trying to do what you did. I got stuck on a bunch of things, most recently this:

my.counts=getBamCounts  (
                bed.frame=hg38exons,
                bam.files=bams,
                include.chr=T
            )

Error in getListElement(x, i, ...) :
GRanges objects don't support [[, as.list(), lapply(), or unlist() at the moment
Calls: getBamCounts ... is.factor -> [[ -> [[ -> getListElement -> getListElement

I believe my GRanges object is well-formed:

head(hg38exons)

GRanges object with 6 ranges and 1 metadata column:
    seqnames        ranges strand |        name
       <Rle>     <IRanges>  <Rle> | <character>
  1     chr1   65419-65433      * |     OR4F5_1
  2     chr1   65520-65573      * |     OR4F5_2
  3     chr1   69037-71585      * |     OR4F5_3
  4     chr1 450740-451678      * |    OR4F29_1
  5     chr1 685716-686654      * |    OR4F16_1
  6     chr1 923923-924948      * |    SAMD11_1
  -------
  seqinfo: 22 sequences from an unspecified genome; no seqlengths

Do you have any idea how I might proceed from here?

Session info:

R version 4.2.2 (2022-10-31)
Platform: x86_64-suse-linux-gnu (64-bit)
Running under: openSUSE Leap 15.5

Matrix products: default
BLAS:   /usr/lib64/R/lib/libRblas.so
LAPACK: /usr/lib64/R/lib/libRlapack.so

locale:
 [1] LC_CTYPE=en_GB.utf8       LC_NUMERIC=C             
 [3] LC_TIME=en_GB.utf8        LC_COLLATE=en_GB.utf8    
 [5] LC_MONETARY=en_GB.utf8    LC_MESSAGES=en_GB.utf8   
 [7] LC_PAPER=en_GB.utf8       LC_NAME=C                
 [9] LC_ADDRESS=C              LC_TELEPHONE=C           
[11] LC_MEASUREMENT=en_GB.utf8 LC_IDENTIFICATION=C      

attached base packages:
[1] stats4    stats     graphics  grDevices utils     datasets  methods  
[8] base     

other attached packages:
[1] ExomeDepth_1.1.16    GenomicRanges_1.50.2 GenomeInfoDb_1.34.9 
[4] IRanges_2.32.0       S4Vectors_0.36.2     BiocGenerics_0.44.0 

loaded via a namespace (and not attached):
 [1] XVector_0.38.0              magrittr_2.0.3             
 [3] zlibbioc_1.44.0             GenomicAlignments_1.34.1   
 [5] BiocParallel_1.32.6         lattice_0.20-45            
 [7] tools_4.2.2                 SummarizedExperiment_1.28.0
 [9] parallel_4.2.2              grid_4.2.2                 
[11] Biobase_2.58.0              matrixStats_1.0.0          
[13] crayon_1.5.2                Matrix_1.6-1               
[15] GenomeInfoDbData_1.2.9      bitops_1.0-7               
[17] codetools_0.2-18            aod_1.3.2                  
[19] RCurl_1.98-1.12             DelayedArray_0.24.0        
[21] compiler_4.2.2              MatrixGenerics_1.10.0      
[23] Rsamtools_2.14.0            Biostrings_2.66.0

[Joel Wallenius]

So it turns out I was just being real dumb. The exon data should be a regular dataframe, not a GRanges object...