SEQanswers

Go Back   SEQanswers > Bioinformatics > Bioinformatics



Similar Threads
Thread Thread Starter Forum Replies Last Post
Bowtie error message nncy1105 Bioinformatics 4 05-24-2013 01:44 PM
HMMSplicer Error Message mjtorresphd Bioinformatics 0 12-13-2012 04:18 PM
Error message in cuffcompare Line RNA Sequencing 3 09-20-2012 12:05 PM
CNVnator : Error message Please Help! KYR Bioinformatics 0 05-20-2012 09:14 PM
GEOquery error message akolman Bioinformatics 8 06-06-2011 07:26 AM

Reply
 
Thread Tools
Old 08-06-2013, 01:10 PM   #1
ToddB
Junior Member
 
Location: Alabama

Join Date: Aug 2013
Posts: 2
Default Error Message in nbinomLRT in DESeq2

Hello all. I am working on some code for my lab and I have run into an error message I can seem to get around. My goal is to test each variable in the model given all other variables are present. Here is the code and the error that are tripping me up.

DESrun=nbinomLRT(DESrun, useOptim=FALSE, maxit = 500 reduced=as.formula(paste("~",paste(Vars[-i],collapse="+"))))

Error in optim(betaRow, objectiveFn, method = "L-BFGS-B", lower = -large, :
L-BFGS-B needs finite values of 'fn'


Vars contain my variables of interest, for example Vars = c("Age","Gender","Disorder"). And DESrun has my full model, dispersion estimates, size factor estimates, and data. I am using DESeq2 version 1.0.18.

If I am still getting this error with the useOptim set to false and maxit set to 500 anyone have a suggestion?

Thanks,
Todd
ToddB is offline   Reply With Quote
Old 08-11-2013, 05:43 AM   #2
Yohann
Junior Member
 
Location: Montréal

Join Date: Aug 2013
Posts: 7
Default

Hi,

I have a pretty similar problem :

Error in optim(betaRow, objectiveFn, method = "L-BFGS-B", lower = -large, :
non-finite finite-difference value [3]


What is weird is that i'm pretty sure it worked with an older version of DESeq2 as i'm re-running an old script.

Did you find a solution ?
Yohann is offline   Reply With Quote
Old 08-11-2013, 07:35 AM   #3
dpryan
Devon Ryan
 
Location: Freiburg, Germany

Join Date: Jul 2011
Posts: 3,478
Default

Can you post the relevant portion of your script (including the design data.frame)? I suspect that something is off in your fit function.
dpryan is offline   Reply With Quote
Old 08-12-2013, 06:34 AM   #4
ToddB
Junior Member
 
Location: Alabama

Join Date: Aug 2013
Posts: 2
Default

No solution yet. But, I have been in contact with Mike Love and he is tracking it down. He sent me the following to see which rows are causing an issue.


if you do:

debug(DESeq2:::fitNbinomGLMs)

this will show which rows are crashing.

nbinomLRT first runs fitNbinomGLMs on the full formula which is not
throwing the error. So you can hit enter through the function once.
The second time it is fitting the reduced formula. On the second time,
after this code chunk, you can identify the rows causing the problem:

# switch based on whether we should also use optim
# on rows which did not converge
if (useOptim) {
rowsForOptim <- which(!betaConv | !rowStable | !rowVarPositive)
} else {
rowsForOptim <- which(!rowStable | !rowVarPositive)
}

The rowsForOptim is the subset of rows which are causing problems.

My current plan is to use this to find the rows causing a problem and remove them, then rerun DESeq2.
ToddB is offline   Reply With Quote
Old 08-14-2013, 04:40 AM   #5
Michael Love
Senior Member
 
Location: Boston

Join Date: Jul 2013
Posts: 333
Default

I now am committing DESeq2 version 1.0.19 to Bioc which uses Nelder-Mead for this back-up optimization rather than L-BFGS-B. This appears to avoid this halting error.

Furthermore, for Todd's data, it is necessary, for now, to also standardize the numeric predictors, e.g.:

# scale the numeric predictors
for (var in c("age","weight")) {
colData(dds)[[var]] <- as.numeric(scale(colData(dds)[[var]]))
}

Then in the end, the log2 fold changes and lfcSE can be multiplied by sd(age) to obtain what would have been the original log2 fold change expected per year for the age variable. (Or you can leave them as they are and they are interpretable as log2 fold changes expected per standard deviation)

I am working on improving the fitting algorithm in the devel branch, so that this manual scaling will not be necessary.
Michael Love is offline   Reply With Quote
Old 08-14-2013, 04:52 AM   #6
Yohann
Junior Member
 
Location: Montréal

Join Date: Aug 2013
Posts: 7
Default

Hi,

I was using a continuous variable that was the number of bacteria as a covariate.
This value scales from ~ 1e4 to 1e7
If i use this scale, i have the error message : "Error in optim(betaRow, objectiveFn, method = "L-BFGS-B", lower = -large, : non-finite finite-difference value [3]"
If i instead divide my values by 1e6 (it's now millions of bacteria), it works.
Yohann is offline   Reply With Quote
Old 08-14-2013, 05:10 AM   #7
Simon Anders
Senior Member
 
Location: Heidelberg, Germany

Join Date: Feb 2010
Posts: 991
Default

Are you sure that you don't want to rather regress on the logarithm of your predictor, given its huge dynamic range?
Simon Anders is offline   Reply With Quote
Old 08-14-2013, 05:14 AM   #8
Yohann
Junior Member
 
Location: Montréal

Join Date: Aug 2013
Posts: 7
Default

Quote:
Originally Posted by Simon Anders View Post
Are you sure that you don't want to rather regress on the logarithm of your predictor, given its huge dynamic range?
That's a good idea, thanks!
I'll check that !
Yohann is offline   Reply With Quote
Old 09-02-2013, 12:55 AM   #9
Michael Love
Senior Member
 
Location: Boston

Join Date: Jul 2013
Posts: 333
Default

I have implemented a change to the fitting algorithm in DESeq2 version >= 1.1.32 which should eliminate the need to standardize predictors as described in my post above.
Michael Love is offline   Reply With Quote
Old 09-04-2013, 10:15 AM   #10
Yohann
Junior Member
 
Location: Montréal

Join Date: Aug 2013
Posts: 7
Default

Hi Michael,

Here is a part of my code :

Code:
colData <- data.frame(
	AFcomp = sub_individuals$AFcomp,
	bacteria_count = sub_individuals$Contact,
	batch = as.factor(paste0("F", sub_individuals$Flowcell)))

dds <- DESeqDataSetFromMatrix(
	countData = sub_reads,
	colData = colData,
	design = ~ AFcomp + bacteria_count + batch)

dds <- estimateSizeFactors(dds)
dds <- estimateDispersions(dds, maxit = 500, quiet = TRUE)
dds <- nbinomLRT(dds, full = design(dds), maxit = 500, reduced = ~ bacteria_count + batch)
Where :
  • AFcomp is a continuous variable scaling from 0 to 1.
  • bacteria_count is a continuous variable from 1e6 to 8e6.
  • batch is a discrete variable (batch name).

if i use the bacteria_count like that, i get this error :

Code:
Error in solve.default(xtwx + ridge) :
  system is computationally singular: reciprocal condition number = 9.67189e-20
If i divide the bacteria_count by 1e6 (now looking at millions of bacteria), it works.

I'm using the 1.1.33 version from the dev branch of bioconductor.

Do you have an idea?

Thanks!

Last edited by Yohann; 09-04-2013 at 10:18 AM.
Yohann is offline   Reply With Quote
Old 09-04-2013, 01:58 PM   #11
Michael Love
Senior Member
 
Location: Boston

Join Date: Jul 2013
Posts: 333
Default

hi Yohann,

I see the problem and will look into it.
Michael Love is offline   Reply With Quote
Old 09-05-2013, 01:22 AM   #12
Michael Love
Senior Member
 
Location: Boston

Join Date: Jul 2013
Posts: 333
Default

hi Yohann,

The problem was coming from some code which goes over rows which did not converge in the standard GLM fitting steps using R's optimization function. So I had not updated this code to account for very large columns in the design matrix. I have now updated this bit of code in version 1.1.35 and added a unit test to confirm that the coefficients and standard errors are identical with or without this extra optimization loop.

If you could check if this solves your error that would be great.
Michael Love is offline   Reply With Quote
Old 09-05-2013, 06:12 AM   #13
Yohann
Junior Member
 
Location: Montréal

Join Date: Aug 2013
Posts: 7
Default

hi Michael, thanks for your answer.

It looks like i can only retrieve the 1.1.34 version on http://bioconductor.org/packages/dev...ml/DESeq2.html.

Is there an other way to get the 1.1.35 ?

Thanks,
Yohann is offline   Reply With Quote
Old 09-05-2013, 06:22 AM   #14
Michael Love
Senior Member
 
Location: Boston

Join Date: Jul 2013
Posts: 333
Default

It should show up on the Bioc site around this time tomorrow.
Michael Love is offline   Reply With Quote
Reply

Thread Tools

Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off




All times are GMT -8. The time now is 01:01 AM.


Powered by vBulletin® Version 3.8.9
Copyright ©2000 - 2018, vBulletin Solutions, Inc.
Single Sign On provided by vBSSO