thanks bill. that's a neat trick that I haven't seen before. now I see what
you're saying
much more clearly. steven: bill's method should be faster than mine because
it won't
have rejection iterations like my idea will.
On Mon, Oct 21, 2013 at 10:52 AM, William Dunlap <wdun...@tibco.com> wrote:
> Try defining the function
> theta345toSigma <- function(theta) {
> cholSigma <- cbind(c(theta[3], 0), c(theta[4], theta[5]))
> crossprod(cholSigma) # t(cholSigma) %*% cholSigma)
> }
> This creates a positive definite matrix for any theta (and
> any positive definite matrix has such a representation, generally
> more than one). It is like using the square root of a quantity
> in the optimizer when you know the quantity must be non-negative.
>
> Then change your
> sigma <- c(theta[3],theta[5],theta[5],theta[4])
> dim(sigma) <- c(2, 2)
> to
> sigma <- theta345toSigma(theta)
>
> If one of your variances is near 0 the optimizer may run into
> trouble at saddlepoints. Others may be able to help better
> with that issue.
>
> Bill Dunlap
> Spotfire, TIBCO Software
> wdunlap tibco.com
>
>
> > -----Original Message-----
> > From: Steven LeBlanc [mailto:ores...@gmail.com]
> > Sent: Sunday, October 20, 2013 9:35 PM
> > To: William Dunlap
> > Subject: Re: [R] nlminb() - how do I constrain the parameter vector
> properly?
> >
> >
> > On Oct 20, 2013, at 6:41 PM, William Dunlap <wdun...@tibco.com> wrote:
> >
> > > Do you mean that your objective function (given to nlminb)
> parameterized
> > > a positive definite matrix, P, as the elements in its upper (or lower)
> triangle?
> > > If so, you could reparameterize it by the non-zero (upper triangular)
> elements
> > > of the Choleski decomposition, C, of P. Compute P as crossprod(C),
> compute
> > > the initial estimate of C as chol(P).
> > >
> > > Bill Dunlap
> > > Spotfire, TIBCO Software
> > > wdunlap tibco.com
> >
> > Hi Bill,
> >
> > I've clipped out the superfluous code to leave the objective function
> and call to nlminb()
> > below.
> >
> > I'm using the first two parameters to construct the vector of means 'u'
> for a bivariate
> > normal, and the final three parameters to construct the corresponding
> covariance matrix
> > 'sigma'. Both are required by dmvnorm(). In summary, nlminb() required a
> vector of
> > parameters, so I supplied the number of parameters I needed nlminb() to
> optimize and
> > simply built the required formats within the function.
> >
> > It didn't occur to me until I saw the error that nlminb() would have no
> way of knowing the
> > proper boundaries of the parameter space unless there is some way to
> communicate the
> > constraints. nlminb() implements simple box constraints, but I don't see
> a way to
> > communicate "parameters 3, 4, and 5 must satisfy 3*4 - 5^2 > 0.
> >
> > Regarding your suggestion, I don't think I understand. Might you
> elaborate?
> >
> > Thanks & Best Regards,
> > Steven
> >
> > > exact<-function(theta,complete,deleted){
> > >>
> > >> u<-c(theta[1],theta[2])
> > >> sigma<-c(theta[3],theta[5],theta[5],theta[4])
> > >> dim(sigma)<-c(2,2)
> > >> -sum(log(dmvnorm(x=complete,mean=u,sigma=sigma)))-
> > >> sum(log(dnorm(one.only,u[1],sigma[1,1])))-
> > >> sum(log(dnorm(two.only,u[2],sigma[2,2])))
> > >> }
> > >>
> >
> nlminb(start=c(0,0,1,1,0),objective=exact,complete=complete,deleted=deleted,control=l
> > >> ist(trace=1))
>
> ______________________________________________
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>
[[alternative HTML version deleted]]
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
