Entering edit mode
Chang Jun Park
▴
10
@chang-jun-park-5086
Last seen 10.2 years ago
Dear Bioconductor developers,
My name is Chris Park who is working as a bioinformatician in
University
Health Network in Toronto. I'm currently working with Illumina
methylation
data by preprocessing them in R, and discovered a problem in which it
could
potentially be a bug.
The bug happens when I call function for Methylation calls:
lumiMethyStatus
{lumi}, which gives this error message:
Error in if (truncate) { : argument is not interpretable as logical
> traceback()
3: gammaFitEM(x, ...)
2: methylationCall(M[, i], ...)
1: lumiMethyStatus(data.lumiMethy.bc.adj.quantile[, c(1:8)])
As shown on traceback(), the error is in gammaFitEM. After looking at
the
code, I found that if iteration reaches the maxInteration (hard-coded
as
50), it will enter:
1 if (iter == maxIteration) {
2 Mode <- c(s[1] + (k[1] - 1) * theta[1], s[2] - (k[2] -
3 1) * theta[2])
4 f1 <- dgamma(x - s[1], shape = k[1], scale = theta[1])
5 f2 <- dgamma(s[2] - x, shape = k[2], scale = theta[2])
6 if (truncate) {
7 f1 <- f1/(pgamma(Mode[2] - s[1], shape = k[1], scale =
theta[1],
8 lower.tail = TRUE))
9 f2 <- f2/(pgamma(s[2] - Mode[1], shape = k[2], scale =
theta[2],
10 lower.tail = TRUE))
11 z1 <- p[1] * f1/(p[1] * f1 + p[2] * f2)
12 z1[x > Mode[2]] <- 0
13 z1[x < Mode[1]] <- 1
14 f1[x > Mode[2]] <- 0
15 f2[x < Mode[1]] <- 0
16 }
17 else {
18 z1 <- p[1] * f1/(p[1] * f1 + p[2] * f2)
19 }
20 z2 <- 1 - z1
21 }
Where on line 6, it checks for logical argument: truncate, but this
logical
variable is nowhere to be found.
What do you suggest is the best solution is?
Theoretically, what is a reason to decide whether to "truncate" or not
when
we reach maximum iteration of the "while" loop?
I am using the latest R version available with the most recent
package, as
shown in sessionInfo():
> sessionInfo()
R version 2.14.1 (2011-12-22)
Platform: x86_64-pc-mingw32/x64 (64-bit)
locale:
[1] LC_COLLATE=English_Canada.1252 LC_CTYPE=English_Canada.1252
[3] LC_MONETARY=English_Canada.1252 LC_NUMERIC=C
[5] LC_TIME=English_Canada.1252
attached base packages:
[1] stats graphics grDevices utils datasets methods base
other attached packages:
[1] KernSmooth_2.23-7
IlluminaHumanMethylation450k.db_1.4.6
[3] org.Hs.eg.db_2.6.4
RSQLite_0.11.1
[5] DBI_0.2-5
AnnotationDbi_1.16.11
[7] lumi_2.6.0
nleqslv_1.9.2
[9] methylumi_2.0.4
Biobase_2.14.0
loaded via a namespace (and not attached):
[1] affy_1.32.1 affyio_1.22.0 annotate_1.32.1
[4] BiocInstaller_1.2.1 grid_2.14.1 hdrcde_2.15
[7] IRanges_1.12.5 lattice_0.20-0 MASS_7.3-16
[10] Matrix_1.0-3 mgcv_1.7-13 nlme_3.1-103
[13] preprocessCore_1.16.0 tools_2.14.1 xtable_1.6-0
[16] zlibbioc_1.0.0
[[alternative HTML version deleted]]
