Error while preprocessing raw data from methylation array using minfi package
1
0
Entering edit mode
Deepa • 0
@f7f3809a
Last seen 7 months ago
Germany

Hi,

I am trying to run the function "preprocessIllumina" from the minfi package on my RGset data. I have been able to import the targets and the RGset. I have also been using the appropriate manifest file and annotation. The code that I am trying to execute is as follows:

preprocessIllumina(RGset, bg.correct = TRUE, normalize = "controls", reference = 1)

Here's the error that is produced:

 Error in h(simpleError(msg, call)) : 
  error in evaluating the argument 'x' in selecting a method for function 'match': object 'CG.controls' not found
In addition: Warning message:
useNames = NA is deprecated. Instead, specify either useNames = TRUE or useNames = TRUE.

Any suggestions on how to fix this?

Best regards
minfi • 1.9k views
ADD COMMENT •
0
Entering edit mode

You need to provide the output of sessionInfo() and from traceback() right after you get the error.

ADD REPLY •
0
Entering edit mode

Here's the sessionInfo and the traceback.

sessionInfo()
R version 4.2.1 (2022-06-23 ucrt)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows 10 x64 (build 19045)

Matrix products: default

locale:
[1] LC_COLLATE=English_United States.utf8  LC_CTYPE=English_United States.utf8    LC_MONETARY=English_United States.utf8
[4] LC_NUMERIC=C                           LC_TIME=English_United States.utf8    

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

other attached packages:
 [1] methylumi_2.44.0                                 FDb.InfiniumMethylation.hg19_2.2.0              
 [3] org.Hs.eg.db_3.16.0                              TxDb.Hsapiens.UCSC.hg19.knownGene_3.2.2         
 [5] GenomicFeatures_1.50.4                           AnnotationDbi_1.60.2                            
 [7] reshape2_1.4.4                                   scales_1.2.1                                    
 [9] QuasR_1.38.0                                     Rbowtie_1.38.0                                  
[11] writexl_1.4.2                                    rlang_1.1.1                                     
[13] BiocParallel_1.32.6                              sesame_1.9.7                                    
[15] sesameData_1.16.0                                ExperimentHub_2.6.0                             
[17] AnnotationHub_3.6.0                              BiocFileCache_2.6.1                             
[19] dbplyr_2.3.3                                     lubridate_1.9.2                                 
[21] forcats_1.0.0                                    stringr_1.5.0                                   
[23] purrr_1.0.2                                      readr_2.1.4                                     
[25] tidyr_1.3.0                                      tibble_3.2.1                                    
[27] ggplot2_3.4.3                                    tidyverse_2.0.0                                 
[29] dplyr_1.1.2                                      GEOquery_2.66.0                                 
[31] HorvathMammalMethylChip40anno.test.unknown_0.2.2 HorvathMammalMethylChip40manifest_0.2.2         
[33] minfi_1.44.0                                     bumphunter_1.40.0                               
[35] locfit_1.5-9.8                                   iterators_1.0.14                                
[37] foreach_1.5.2                                    Biostrings_2.66.0                               
[39] XVector_0.38.0                                   SummarizedExperiment_1.28.0                     
[41] Biobase_2.58.0                                   MatrixGenerics_1.10.0                           
[43] matrixStats_1.0.0                                GenomicRanges_1.50.2                            
[45] GenomeInfoDb_1.34.9                              IRanges_2.32.0                                  
[47] S4Vectors_0.36.2                                 BiocGenerics_0.44.0                             
[49] RColorBrewer_1.1-3                              

loaded via a namespace (and not attached):
  [1] plyr_1.8.8                    splines_4.2.1                 digest_0.6.33                 htmltools_0.5.6              
  [5] fansi_1.0.4                   magrittr_2.0.3                memoise_2.0.1                 BSgenome_1.66.3              
  [9] tzdb_0.4.0                    limma_3.54.2                  annotate_1.76.0               askpass_1.1                  
 [13] timechange_0.2.0              siggenes_1.72.0               prettyunits_1.1.1             jpeg_0.1-10                  
 [17] colorspace_2.1-0              blob_1.2.4                    rappdirs_0.3.3                crayon_1.5.2                 
 [21] RCurl_1.98-1.12               genefilter_1.80.3             VariantAnnotation_1.44.1      survival_3.3-1               
 [25] glue_1.6.2                    gtable_0.3.4                  zlibbioc_1.44.0               DelayedArray_0.24.0          
 [29] wheatmap_0.2.0                Rhdf5lib_1.20.0               HDF5Array_1.26.0              DBI_1.1.3                    
 [33] rngtools_1.5.2                Rcpp_1.0.11                   xtable_1.8-4                  progress_1.2.2               
 [37] bit_4.0.5                     mclust_6.0.0                  preprocessCore_1.60.2         httr_1.4.7                   
 [41] ellipsis_0.3.2                pkgconfig_2.0.3               reshape_0.8.9                 XML_3.99-0.14                
 [45] deldir_1.0-9                  DNAcopy_1.72.3                utf8_1.2.3                    tidyselect_1.2.0             
 [49] later_1.3.1                   munsell_0.5.0                 BiocVersion_3.16.0            tools_4.2.1                  
 [53] cachem_1.0.8                  cli_3.6.1                     generics_0.1.3                RSQLite_2.3.1                
 [57] fastmap_1.1.1                 yaml_2.3.7                    bit64_4.0.5                   beanplot_1.3.1               
 [61] scrime_1.3.5                  randomForest_4.7-1.1          KEGGREST_1.38.0               nlme_3.1-157                 
 [65] doRNG_1.8.6                   sparseMatrixStats_1.10.0      mime_0.12                     nor1mix_1.3-0                
 [69] xml2_1.3.5                    biomaRt_2.54.1                compiler_4.2.1                rstudioapi_0.15.0            
 [73] filelock_1.0.2                curl_5.0.2                    png_0.1-8                     interactiveDisplayBase_1.36.0
 [77] stringi_1.7.12                GenomicFiles_1.34.0           lattice_0.20-45               Matrix_1.6-1                 
 [81] multtest_2.54.0               vctrs_0.6.3                   pillar_1.9.0                  lifecycle_1.0.3              
 [85] rhdf5filters_1.10.1           BiocManager_1.30.19           data.table_1.14.8             bitops_1.0-7                 
 [89] httpuv_1.6.11                 rtracklayer_1.58.0            latticeExtra_0.6-30           hwriter_1.3.2.1              
 [93] R6_2.5.1                      BiocIO_1.8.0                  ShortRead_1.56.1              promises_1.2.1               
 [97] codetools_0.2-18              MASS_7.3-57                   rhdf5_2.42.1                  openssl_2.1.0                
[101] rjson_0.2.21                  withr_2.5.0                   GenomicAlignments_1.34.1      Rsamtools_2.14.0             
[105] GenomeInfoDbData_1.2.9        hms_1.1.3                     quadprog_1.5-8                grid_4.2.1                   
[109] base64_2.0.1                  DelayedMatrixStats_1.20.0     illuminaio_0.40.0             shiny_1.7.5                  
[113] interp_1.1-4                  restfulr_0.0.15

traceback()
9: h(simpleError(msg, call))
8: .handleSimpleError(function (cond) 
   .Internal(C_tryCatchHelper(addr, 1L, cond)), "object 'CG.controls' not found", 
       base::quote(match(CG.controls, rownames(Green))))
7: match(CG.controls, rownames(Green))
6: matrixStats::colMeans2(x, rows = rows, cols = cols, na.rm = na.rm, 
       dim. = dim., ..., useNames = useNames)
5: .local(x, rows, cols, na.rm, ..., useNames = useNames)
4: colMeans2(Green, rows = match(CG.controls, rownames(Green)))
3: colMeans2(Green, rows = match(CG.controls, rownames(Green)))
2: normalize.illumina.control(rgSet, reference = reference)
1: preprocessIllumina(RGset, bg.correct = TRUE, normalize = "controls", 
       reference = 1)
ADD REPLY •
0
Entering edit mode

What array is this? One of those Horvath mammalian arrays, or an EPIC?

ADD REPLY •
0
Entering edit mode

Could you show the output of RGset ? What I suspect is that you have an EPICv2 array and it seems minfi does not support this type of array for the moment

ADD REPLY •
0
Entering edit mode

The idat files are generated from the HorvathMammalMethylChip40. I provide the correct manifest file in subsequent steps during the analysis.

class: RGChannelSet 
dim: 404407 12 
metadata(0):
assays(2): Green Red
rownames(404407): 1600149 1600245 ... 99810932 99810936
rowData names(0):
colnames(12): 206648980005_R01C01 206648980005_R02C01 ... 206648980005_R05C02 206648980005_R06C02
colData names(15): OriginalOrderInBatch Sample_Name ... Basename filenames
Annotation
  array: Unknown
  annotation: Unknown
ADD REPLY •
0
Entering edit mode

Sure it does. But you need to get the manifest from Zuguang Gu https://jokergoo.github.io/IlluminaHumanMethylationEPICv2manifest/

> library(BiocManager)
> install("jokergoo/IlluminaHumanMethylationEPICv2manifest")
> install("jokergoo/IlluminaHumanMethylationEPICv2anno.20a1.hg38")
> RGset <- read.metharray.exp(targets = targets)
> annotation(RGset <- c(array = "IlluminaHumanMethylationEPICv2", annotation = "20a1.hg38")
ADD REPLY •
0
Entering edit mode
@james-w-macdonald-5106
Last seen 3 hours ago
United States

The Horvath arrays may not have control probes. You can check. Here's what is on an EPICV2

> z <- getProbeInfo(dat, "Control")
> table(z$Type)

 BISULFITE CONVERSION I BISULFITE CONVERSION II               EXTENSION 
                     10                       4                       4 
          HYBRIDIZATION                NEGATIVE         NON-POLYMORPHIC 
                      3                     411                       9 
                 NORM_A                  NORM_C                  NORM_G 
                     27                      58                      27 
                 NORM_T             RESTORATION           SPECIFICITY I 
                     58                       1                      12 
         SPECIFICITY II                STAINING          TARGET REMOVAL 
                      3                       6                       2

And if you don't see NORM_A - NORM_T in that list, then that's your problem.
ADD COMMENT •
0
Entering edit mode

I do see NORM_A NORM_T in the list

z <- getProbeInfo(manifest, "Control")
> table(z$Type)

 BISULFITE CONVERSION I BISULFITE CONVERSION II               EXTENSION           HYBRIDIZATION                NEGATIVE 
                     10                      10                       4                       3                     787 
        NON-POLYMORPHIC                  NORM_A                  NORM_C                  NORM_G                  NORM_T 
                      4                      48                      58                      44                      61 
          SPECIFICITY I          SPECIFICITY II                STAINING          TARGET REMOVAL 
                     12                       6                       4                       2
ADD REPLY •
0
Entering edit mode

It appears that you are getting that from the manifest object. What do you get when you do

z <- getProbeInfo(RGset, "Control")
ADD REPLY •
0
Entering edit mode

You might need to do

annotation(RGset) <- c(array = "HorvathMammalMethylChip40manifest", annotation = "test.unknown")
ADD REPLY •
0
Entering edit mode

Yes, I do this in the subsequent steps

ADD REPLY •
0
Entering edit mode

Here's what it looks like

y <- getProbeInfo(RGset1, "Control")
> table(y$Type)

 BISULFITE CONVERSION I BISULFITE CONVERSION II               EXTENSION           HYBRIDIZATION                NEGATIVE 
                     10                       1                       4                       1                     396 
        NON-POLYMORPHIC                  NORM_A                  NORM_C                  NORM_G                  NORM_T 
                      1                      27                      58                      27                      58 
          SPECIFICITY I          SPECIFICITY II                STAINING          TARGET REMOVAL 
                      4                       3                       4                       2
ADD REPLY •
0
Entering edit mode

Oh wait. The code for normalize.illumina.control has this

if (.is450k(rgSet) || .isEPIC(rgSet) || .isAllergy(rgSet)) {
        AT.controls <- getControlAddress(object = rgSet, controlType = c("NORM_A", 
            "NORM_T"))
        CG.controls <- getControlAddress(object = rgSet, controlType = c("NORM_C", 
            "NORM_G"))
    }
    if (.is27k(rgSet)) {
        AT.controls <- getControlAddress(object = rgSet, controlType = "Normalization-Red")
        CG.controls <- getControlAddress(object = rgSet, controlType = "Normalization-Green")
    }
    Green.avg <- colMeans2(Green, rows = match(CG.controls, rownames(Green)))
    Red.avg <- colMeans2(Red, rows = match(AT.controls, rownames(Red)))

And your array is going to fail all those .isXXX tests, which are going to do tests like

annotation(object)["array"] == "IlluminaHumanMethylationEPIC"

So maybe you need to first change the array to e.g., IlluminaHumanMethylationEPIC and change back later. No guarantee that will work though. You might contact the Horvath group or look through their GitHub to see if there is a workaround that they recommend.

ADD REPLY •
0
Entering edit mode

Thank you very much! I shall try both to check what works.

ADD REPLY •

Login before adding your answer.

Similar Posts
Loading Similar Posts
Traffic: 583 users visited in the last hour
Content Search
Users
Tags
Badges
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6