Re: [R] question about nlminb

From: Spencer Graves <spencer.graves_at_pdf.com>
Date: Wed, 09 Apr 2008 18:30:40 -0700

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.
>
>



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. Received on Thu 10 Apr 2008 - 01:46:54 GMT

Archive maintained by Robert King, hosted by the discipline of statistics at the University of Newcastle, Australia.
Archive generated by hypermail 2.2.0, at Sun 13 Apr 2008 - 21:30:29 GMT.

Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-help. Please read the posting guide before posting to the list.

list of date sections of archive