 Email discussion lists for the UK Education and Research communities  ## BUGS@JISCMAIL.AC.UK

#### View:

 Message: [ First | Previous | Next | Last ] By Topic: [ First | Previous | Next | Last ] By Author: [ First | Previous | Next | Last ] Font: Proportional Font  LISTSERV Archives  BUGS Home  BUGS 2006

#### Options  Subscribe or Unsubscribe   Log In   Get Password Subject: Non-convergence in nonlinear parameter estimation problem

From:  Date: Wed, 2 Aug 2006 10:35:30 -0500

Content-Type: text/plain

Parts/Attachments:  text/plain (106 lines)
 ```Hello all, I am trying to use WinBUGS to estimate the parameters of a nonlinear function, from 'data' with known Gaussian errors. (These 'data' are simulated). The function is of the form: f(x) = aa + bb*sin(cc*x^2 + dd*x + ee)   aa = 2.0   bb = 3.0   cc = 1.5   dd = -4.0   ee = 3.0 My 'data' consists of 30 (x,y) pairs. To generate each pair, I use the following procedure:    1. choose x[i] from a uniform distribution on the range (0,5)    2. set y[i] = f(x[i]) + e[i] where e[i] is selected from a Gaussian distribution with mean 0 and standard deviation 0.2 I then use WinBUGS, through R, to estimate the parameters of the function f(x) from the data. My BUGS model is expressed in the file model_1.txt: ----------------%<-------------------------------------------------- # File model_1.txt # var x1[NDATA], y1[NDATA], yfit[NDATA], aa, bb, cc, dd, ee, tau, sigma model {    for (i in 1:NDATA)    {      yfit[i] <- aa + bb * sin(cc*x1[i]*x1[i] +dd*x1[i] + ee);      y1[i] ~ dnorm(yfit[i], tau);    }    tau <- 25.0    aa ~ dnorm(2.0, 0.01);    bb ~ dnorm(3.0, 0.01);    cc ~ dnorm(1.5, 0.01);    dd ~ dnorm(-4.0, 0.01);    ee ~ dnorm(3.0, 0.01); } ----------------%<-------------------------------------------------- I run the calculation from R, using the R script model_1.R, below. The variables NDATA, x1 and y1 contain the number of points (30), and the x and y values of my simulated data. ----------------%<-------------------------------------------------- # File model_1. # data <- list("NDATA", "x1", "y1") inits <- function()           {             aa = rnorm(1, 2.0, 5.0)             bb = rnorm(1, 3.0, 5.0) cc = rnorm(1, 1.5, 5.0) dd = rnorm(1, -4.0, 5.0) ee = rnorm(1, 3.0, 5.0)           } parameters = c("aa", "bb", "cc", "dd", "ee") model.file = "model_1.txt" fit1.sim <- bugs( data                  , inits=NULL                  , parameters                  , model.file                  , n.chains=5                  , n.iter=100*1000                  , coda=TRUE                  , debug=TRUE                  , DIC=TRUE                  ) cda1 <- read.bugs(fit1.sim) ----------------%<-------------------------------------------------- Even with what I think is a very large number of iterations (100,000 -- which takes about 30 minutes on my laptop), it seems that the chains have not converged. The resulting posterior distributions for each of the five parameters is multimodal, and the trace plots (produced by plot(cda1) in R) show no hint of convergence. I do not understand what I might be doing wrong in my use of WinBUGS, and would appreciate any assistance someone with more expertise might provide. best regards, Dr. Marc Paterno Fermi National Accelerator Laboratory ------------------------------------------------------------------- This list is for discussion of modelling issues and the BUGS software. For help with crashes and error messages, first mail [log in to unmask] To mail the BUGS list, mail to [log in to unmask] Before mailing, please check the archive at www.jiscmail.ac.uk/lists/bugs.html Please do not mail attachments to the list. To leave the BUGS list, send LEAVE BUGS to [log in to unmask] If this fails, mail [log in to unmask], NOT the whole list ```