[R] logistic transformation using nlminb
Dear all,
I want to find the optimal values of a vector, x (with 6 elements)
say, satisfying the following conditions:
1. for all x=0
2. sum(x)=1
3. x[5]=0.5 and x[6]=0.5
For the minimisation I'm using nlminb and to satisfy the first 2
conditions the logistic transformation is used with box constraints
for condition 3. However, I don't seem to be able to get the values x
that i'm expecting from fitting a simpler model. For example,
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(1, mean=c(3, 2.75, 2.8, 3, 2.9), sd=c(0.1,
0.1,0.12, 0.15, 0.10), cov=corr)
# here is the simple optimisation function that allows the 5th
element to be potentially negative
obj.funA - function(opt, mat) {
opt - c(opt, 1-sum(opt))
LinearComb - mat%*%opt
obj - -min(LinearComb)
obj
}
# I want to put an upper boundary constraint on the first variable -
not being greater than 0.35 and the rest being between 0 and 1
opt - nlminb(rep(0,4), lower=rep(0,4), upper=c(0.35, 1, 1, 1) ,
obj.funA, mat=mat)
opt.x - opt$parameters
opt.x - c(opt.x, 1-sum(opt.x))
opt.x
# Result
[1] 0.3402 0.06475651 0. 0.16561760 0.41962686
However, I don't get the same result from the logistic transformation
obj.funB - function(opt, mat) {
opt - c(opt, 0)
opt - exp(opt)/sum(exp(opt))
LinearComb - mat%*%opt
obj - -min(LinearComb)
obj
}
opt - nlminb(rep(0,4), upper=c(log(0.35), NA, NA, NA), obj.funB, mat=mat)
opt.x - opt$parameters
opt.x - c(opt.x, 0)
opt.x - exp(opt.x)/sum(exp(opt.x))
opt.x
# Result
[1] 2.355325e-001 6.339398e-009 1.202751e-004 9.139718e-002 6.729500e-001
I don't know how to obtain the same answer for both optimisations. In
reality, my own optimisation typically gives negative values for the standard
minimisation- so I have no choice but to use a more advanced method.
Also, there appears to be a dependency between the first and the last,
i.e. 2.355325e-001/6.729500e-001 =0.35
Does anyone know why the logistic doesn't give the same answer as the
simpler method?
Any help is much appreciated.
Regards,
John
__
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.
[R] logistic transformation using nlminb
Dear all,
I want to find the optimal values of a vector, x (with 6 elements)
say, satisfying the following conditions:
1. for all x=0
2. sum(x)=1
3. x[5]=0.5 and x[6]=0.5
For the minimisation I'm using nlminb and to satisfy the first 2
conditions the logistic transformation is used with box constraints
for condition 3. However, I don't seem to be able to get the values x
that i'm expecting from fitting a simpler model. For example,
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(1, mean=c(3, 2.75, 2.8, 3, 2.9), sd=c(0.1,
0.1,0.12, 0.15, 0.10), cov=corr)
# here is the simple optimisation function that allows the 5th
element to be potentially negative
obj.funA - function(opt, mat) {
opt - c(opt, 1-sum(opt))
LinearComb - mat%*%opt
obj - -min(LinearComb)
obj
}
# I want to put an upper boundary constraint on the first variable -
not being greater than 0.35 and the rest being between 0 and 1
opt - nlminb(rep(0,4), lower=rep(0,4), upper=c(0.35, 1, 1, 1) ,
obj.funA, mat=mat)
opt.x - opt$parameters
opt.x - c(opt.x, 1-sum(opt.x))
opt.x
# Result
[1] 0.3402 0.06475651 0. 0.16561760 0.41962686
However, I don't get the same result from the logistic transformation
obj.funB - function(opt, mat) {
opt - c(opt, 0)
opt - exp(opt)/sum(exp(opt))
LinearComb - mat%*%opt
obj - -min(LinearComb)
obj
}
opt - nlminb(rep(0,4), upper=c(log(0.35), NA, NA, NA), obj.funB, mat=mat)
opt.x - opt$parameters
opt.x - c(opt.x, 0)
opt.x - exp(opt.x)/sum(exp(opt.x))
opt.x
# Result
[1] 2.355325e-001 6.339398e-009 1.202751e-004 9.139718e-002 6.729500e-001
I don't know how to obtain the same answer for both optimisations. In
reality, my own optimisation typically gives negative values for the standard
minimisation- so I have no choice but to use a more advanced method.
Also, there appears to be a dependency between the first and the last,
i.e. 2.355325e-001/6.729500e-001 =0.35
Does anyone know why the logistic doesn't give the same answer as the
simpler method?
Any help is much appreciated.
Regards,
John
__
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.
Re: [R] question about nlminb
Hi Spencer,
Thanks for your email.
Do you have a reference for generating the variance-covariance matrix
from the restricted variance-covariance? Is this a well known
technique?
Regards,
John
On 10/04/2008, Spencer Graves [EMAIL PROTECTED] wrote:
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(1, 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(1, 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 [EMAIL PROTECTED] 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@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@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@r-project.org mailing
Re: [R] question about nlminb
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(1, 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 [EMAIL PROTECTED] 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@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@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] question about nlminb
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@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] question about nlminb
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@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.
