On Sunday, November 3, 2013 3:17 PM, IZHAK shabsogh wrote:
Hi,
another problem after getting the above hessian matrix i am suppose to use it and fine the sum of the product
of the hessian and the residuals, for the hessian is generated from the above while the residual is also obtain from
the code below :
freg<-function(){
+ reg<- lm(y ~ x1 + x2 , data=x)
+ fit<-nls(y ~ x1 / (1+ b1*x2^b2),data = xx, start =
+ list(b1= 0.0215,b2=-0.1675 ))
+ residuals<-data.frame(residuals(fit))
+ rresult<-list(reg,fit,residuals)
+ print(rresult)
+ }
> freg()
for me to do that i try this code to get the sum of the product of the residuals and the hessian as
residual = e, hessian =h.
Matrix = e1*h1 +
e2*h2+ e3*h3 + . . .e29h29
he <- function(){
for(i in 1:29){
tkp<-hessianList[[i]]
tkp1<-residuals[i]
kpo<-tkp*tkp1
print(kpo)
}
}
he()
my problem is on how to evaluate the above that is (Matrix = e1*h1 + e2*h2+ e3*h3 + . . .e29h29)
kindly guide me on how to go about this problem also if there is any mistake in the flow of my code
need your correction please because i am a beginner
thanks your understanding
On Saturday, November 2, 2013 1:24 AM, Dennis Murphy wrote:
Hi:
In R, it is common to write a function to do the hard work for a
generic single instance and then use an 'apply family' function to map
the function to a set of instances. In this case, the mapply()
function fits your problem since you have multiple input vectors of
the same length.
# Function to generate a 2 x 2 matrix taking x1 and x2 as
# scalar input arguments, using your given values of b1 and b2
# as default values. They can be changed when the function
# is called.
hessianGen <- function(x1, x2, b1 = 4.286, b2 = 1.362)
{
gh<-matrix(0,2,2)
exp0<-(1+b1*x2^b2)
exp1<-x1*x2^b2*log(x2)
exp3<-x1*b1*x2^b2*(log(x2))^2
gh[1,1]<-2*x2^(2*b2)*exp0/exp0^4
gh[1,2]<--(exp0^2*exp1 - 2*b1*x2^b2*exp0*exp1)/exp0^4
gh[2,1]<--(exp3*exp0^2-2*exp0*b1^2*x2^b2*log(x2)*exp1)/exp0^4
gh[2,2]<--(exp1*exp0^2-2*exp0*x2^b2*b1*exp1)/exp0^4
gh
}
# Use mapply to recursively run hessianGen() on the corresponding
# elements of x1 and x2. The argument SIMPLIFY = FALSE is
# required to ensure that the matrix output of the function is kept.
# Returns a list of matrices.
hessianList <- mapply(hessianGen, x1, x2, SIMPLIFY = FALSE)
> length(hessianList)
[1] 29
> hessianList[[1]]
[,1] [,2]
[1,] 0.004186016 0.2470583
[2,] -2.104838360 0.2470583
> hessianList[[2]]
[,1] [,2]
[1,] 0.006006215 0.05384947
[2,] 0.211386483 0.05384947
Your definition of gh[1, 2] and gh[2, 2] are the same. If that's the
intent, you could save some computation in the function and just write
gh[2, 2] <- gh[1, 2]; if not, you might want to look more closely at
your definitions.
Dennis
On Fri, Nov 1, 2013 at 3:06 AM, IZHAK shabsogh wrote:
> below is a code to compute hessian matrix , which i need to generate 29 number of different matrices for example first element in x1 and x2 is use to generate let say matrix (M1) and second element in x1 and x2 give matrix (M2) upto matrix (M29) corresponding to the total number of observations and b1 and b2 are constant.
> can some one guide me or help to implement this please. I did not understand how to construct the loop which i think it should be
>
> about 3 dfferent loops
> example i = 1 to 29 number of matrices
> j1=1 t0 2 row of matirx
> j2= 1 to 2 ncol of matrix
>
>
> x1<-c(5.548,4.896,1.964,3.586,3.824,3.111,3.607,3.557,2.989,18.053,3.773,1.253,2.094,2.726,1.758,5.011,2.455,0.913,0.890,2.468,4.168,4.810,34.319,1.531,1.481,2.239,4.204,3.463,1.727)
> y<-c(2.590,3.770,1.270,1.445,3.290,0.930,1.600,1.250,3.450,1.096,1.745,1.060,0.890,2.755,1.515,4.770,2.220,0.590,0.530,1.910,4.010,1.745,1.965,2.555,0.770,0.720,1.730,2.860,0.760)
> x2<-c(0.137,2.499,0.419,1.699,0.605,0.677,0.159,1.699,0.340,2.899,0.082,0.425,0.444,0.225,0.241,0.099,0.644,0.266,0.351,0.027,0.030,3.400,1.499,0.351,0.082,0.518,0.471,0.036,0.721)
> b1<-4.286b2<-1.362
>
> n<-29
> for(i in 1:n){
> gh<-matrix(0,2,2)
>
exp0<-(1+b1*x2^b2)
> exp1<-x1*x2^b2*log(x2)
> exp3<-x1*b1*x2^b2*(log(x2))^2
> gh[1,1]<-2*x2^(2*b2)*exp0/exp0^4
> gh[1,2]<--(exp0^2*exp1 - 2*b1*x2^b2*exp0*exp1)/exp0^4
> gh[2,1]<--(exp3*exp0^2-2*exp0*b1^2*x2^b2*log(x2)*exp1)/exp0^4
> gh[2,2]<--(exp1*exp0^2-2*exp0*x2^b2*b1*exp1)/exp0^4
> }
> [[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.
>
[[alternative HTML version deleted]]