CODEWITHSUNDEEP
rmatrixregression
hyat
hyat

Reputation: 1057

how to perform nls fro three matrix?

If we have these matrix as example:

    x=matrix(runif(50 ),5,5)
    z=matrix(runif(60 ),5,5)
    y=exp(0.5*x+0.2*z+0.3)^1

the values in these three matrix change form one day to another so I have two years (daily,x1,x2,x3,etc)of these matrix(speaking of my real data) to fit this equation y=exp(s1x+s2x+s3)^s4 use normal Nls as: 

        fit <- nls(y ~ exp(s1*x+s2*z+s3)^s4,
             data = data,
             start = list(s1 = 0.5, s2 = 0.2, s3 = 0.3,s4=1))

but my problem is that I want to do this fit for all corresponding values in three matrix for example start taking values from:

                   > use all values of x1[1,1], x2[1,1],.....etc
                  [1] 0.3617776 .......etc
                   > and the corresponding form z1[1,1], z2[1,1].....etc
                 [1] 0.5544851 .......etc
                   > y1[1,1], y2[1,1],etc .....
                 [1] 1.807213 .......etc

Find the fit parameters then do the same for all other corresponding values. Finally I get a matrix with, for example, values of s1, another matrix for s2 etc….

Any help is appreciated(Note: I have to use it as it is and not to linearise it or use lm.fit or lm)

Upvotes: 1

Views: 110

Answers (1)

pbible
pbible

Reputation: 1259

I think the x,y,z as matrices are throwing you off. As stated by MrFlick in the comment you need multiple observations to perform the regression. Here is an example of nlm applied to a simple 3d normal surface with some noise.

You need a n-by-3 matrix of observations in the 3d space. Here I make some noisy points along a Gaussian curve and rotate it to form a surface.

This example uses the rgl package.

library(rgl)

z_rotate <- function(mat,rads){
    rot.mat <- matrix(c(cos(rads),-sin(rads),0,
                    sin(rads),cos(rads),0,
                    0, 0, 1),nrow=3,ncol=3)
    mat %*% rot.mat
}

x <- seq(-2,2,0.1)
x <- x + rnorm(length(x),sd=0.5)
y <- seq(-2,2,0.1)
#y <- y + rnorm(length(y),sd=0.5)
z <-  0.5*exp(-x^2/2 + (-y^2/2)) + rnorm(1,sd=0.6)

m <- as.matrix(cbind(x,y,z))

now the matrix m is just the initial points about the normal curve.

points <- Reduce(rbind,lapply(1:8,function(n){z_rotate(m,n*(pi/8))}))
colnames(points) <- c("x","y","z")

The previous command just calls the rotate function then appends the results as new rows in the matrix.

Now the points matrix represents the observations to be approximated with our nlm function.

> dim(points)
[1] 328   3

Now we can view these points in 3d with RGL.

plot3d(points[,1],points[,2],points[,3],type='s',size=0.3)

enter image description here

Use nlm to fit a function that approximates this data.

fit <- nls(z ~ s1*exp(-s2*x^2 + (-s3*y^2)) +s4,
        data = data.frame(points),
        start = list(s1 = 0.3, s2 = 0.3, s3 = 0.3, s4=0.3),
        #trace=TRUE,
        control = list(warnOnly=TRUE))

It produces this fit:

> fit
Nonlinear regression model
  model: z ~ s1 * exp(-s2 * x^2 + (-s3 * y^2)) + s4
   data: data.frame(points)
     s1      s2      s3      s4 
 0.5000  0.5000  0.5000 -0.1084 
 residual sum-of-squares: 1.493e-31

Number of iterations till stop: 50 
Achieved convergence tolerance: 0.05046
Reason stopped: number of iterations exceeded maximum of 50

We can see that it recovers the values we had, namely 0.5 for the coefficients when we compare the formulas.

z <-  0.5*exp(-x^2/2 + (-y^2/2)) + rnorm(1,sd=0.6)

z ~ s1*exp(-s2*x^2 + (-s3*y^2)) +s4

I know its a toy example but it might give you some ideas as to how to set up your data. I think you will need to structure your data in an n-by-3 observations matrix to be able to do your regression with nlm.

Upvotes: 2

Related Questions