[R] question about nlminb

[Spencer Graves]

Hi, John: 

      I just got the following error right after the attempt to use 
'rmvnorm'. 

Error: could not find function "rmvnorm"

      I tried 'library(mvtnorm)', but the 'rmvnorm' in that package gave 
me the following: 

Error in rmvnorm(10000, mean = c(3, -20, -10, 3, 2), sd = c(0.1, 15, 4,  :
  unused argument(s) (sd = c(0.1, 15, 4, 0.15, 0.8), cov = c(1, 0.5, 
0.5, 0.5, 0.5, 0.5, 1, 0.5, 0.5, 0.5, 0.5, 0.5, 1, 0.5, 0.5, 0.5, 0.5, 
0.5, 1, 0.5, 0.5, 0.5, 0.5, 0.5, 1))

      Then I got 87 hits to RSiteSearch("rmvnorm", 'fun'). 

      Regarding, "How would I go about getting the entire 
variance-covariance matrix?", this is really easy:  Just define a 5 x 4 
matrix A so the 5-dimensional 'opt' you want is a constant plus A times 
the 4-dimensional restricted 'opt'.  [Please don't use the same name for 
two different things like 'opt' in your function below.  Half a century 
ago, when Fortran was young, that was very smart programming.  Today, 
it's primary effect it to make it difficult for mere mortals to 
understand your code.  Use something like 'opt5' for one and 'opt4' for 
the other.] 

      Then

      Sig5 = A %*% vcov.nlminb(opt) %*% t(A). 

      I believe the A matrix you want is as follows: 

A = matrix(NA, dim=c(5, 4))
A[1:4, 1:4] <- diag(4)
A[5, ] <- (-1)

      However, since your example was not self contained, I have not 
tested this.     
    
      Hope this helps. 
      Spencer

John Pitchard wrote:
> Hi Spencer,
>
> Sorry for not producing code as a worked example.
>
>
> Here's an example:
> ==================================
> # setting the seed number
> set.seed(0)
> # creating a correlation matrix
> corr <- diag(5)
> corr[lower.tri(corr)] <- 0.5
> corr[upper.tri(corr)] <- 0.5
>
> # Data for the minimisation
> mat <- rmvnorm(10000, mean=c(3, -20, -10, 3, 2), sd=c(0.1, 15, 4,
> 0.15, 0.8), cov=corr)
>
> obj.fun <- function(opt, mat) {
>    opt <- c(opt, 1-sum(opt))
>    LinearComb <- mat%*%opt
>    obj <- -min(LinearComb)
>    obj
> }
>
> opt <- nlminb(rep(0,4), lower=rep(-3, 4), upper=rep(3, 4), obj.fun, mat=mat)
> opt.x <- opt$parameters
> opt.x <- c(opt.x, 1-sum(opt.x))
>
> # using vcov.nlminb from the MASS library to obtain the covariance matrix
> vcov.nlminb(opt)
> ====================================
>
>
> I have a variance-covariance matrix for 4 of the elements in the
> vector but not the last component. How would I go about getting the
> entire variance-covariance matrix?
>
> Thanks in advance for any help.
>
> Regards,
> John
>
>
>
>
> On 09/04/2008, Spencer Graves <spencer.graves at pdf.com> wrote:
>   
>>     Have you considered optimizing over x1 = x[1:(length(x)-1]?   You could feed a wrapper function 'f2(x1, ...)' that computes xFull = c(x1, 1-sum(x1)) and feeds that to your 'fn'.
>>     If this makes sense, great.  Else, if my answer is not useful, be so kind as to PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
>>     Spencer
>>
>> John Pitchard wrote:
>>
>>     
>>>  Dear All,
>>>
>>> I wanted to post some more details about the query I sent to s-news last
>>> week.
>>>
>>> I have a vector with a constraint. The constraint is that the sum of the
>>> vector must add up to 1 - but not necessarily positive, i.e.
>>>
>>> x[n] <- 1 -(x[1] + ...+x[n-1])
>>>
>>> I perform the optimisation on the vector x such that
>>>
>>> x <- c(x, 1-sum(x))
>>>
>>> In other words,
>>>
>>> fn <- function(x){
>>>  x <- c(x, 1 - sum(x))
>>>  # other calculations here
>>> }
>>>
>>> then feed this into nlminb()
>>>
>>> out <- nlminb(x, fn)
>>> out.x <- out$parameters
>>> out.x <- c(out.x, 1 - sum(out.x))
>>> out.x
>>>
>>> I would like to calculate standard errors for each of the components of x.
>>> Is this possible by outputing the Hessian matrix? Furthermore, how would I
>>> calculate this for the last component (if this is indeed possible) which has
>>> the restriction (i.e. 1-sum(out.x))?
>>>
>>> Any help would be much appreciated.
>>>
>>> Regards,
>>> John
>>>
>>>        [[alternative HTML version deleted]]
>>>
>>> ______________________________________________
>>> R-help at 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.
>>>
>>>
>>>       
>
> ______________________________________________
> R-help at 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.
>
>