CODEWITHSUNDEEP
r
THIOYE
THIOYE

Reputation: 11

Non linear regression with nls2 package

I work on agglomeration effect and i want to run a non linear regression with package nls2.

I am trying to run this model with R

dens=runif(100)
surf=rnorm(100, 10, 2)  
zone=seq(1,100,1)  
donnees<-data.frame(dens,surf,zone)   
attach(donnees) 
donnees$salaire<-rnorm(100, 1000,3)    
mp<-rep(0,100)    
MP<-rep(0,100)
MPfonc<-function(alpha){    
for (i in 1:100){    
   for (j in 1:100){    
     if(j!=i){    
    mp[j]<- dens[j]/(surf[i]-surf[j])^alpha
     }
  }    
  MP[i]<-sum(mp)    
}    
  return(MP)    
}

fo <- salaire ~ const+ gamma1*dens+gamma2*surf+gamma3*MPfonc(alpha)

gstart <- data.frame(const = c(-100, 100), gamma1 = c(-10, 10),
       gamma2 = c(-10, 10),gamma3 = c(-10, 10), alpha=c(-10, 10))
fm <- nls2(fo, start = gstart,  alg = "plinear-random")

It does not run and I think it is a problem of alpha. Can nls2 function accept a function (MP(alpha)) as an input?

Here is the specification of my model:

nonlinear regression equations

Upvotes: 0

Views: 1007

Answers (1)

G. Grothendieck
G. Grothendieck

Reputation: 269371

The problems are:

  • set.seed should be set to make the code reproducible

  • salaire is not defined -- it is defined in the data frame donnees but donnees is never used after that.

  • the elements summed in the sum call in MPfonc include NaN or NA elements so the sum becomes similarly undefined

  • For plinear algorithms the RHS of the formula must evaluate to a matrix of coefficients of the linear parameters.

  • For the plinear algorithms provide starting values only for the non-linear parameters (i.e. only for alpha).

  • the nls2 package is never loaded. A library statement is needed.

  • code posted to SO should be indented 4 spaces to cause it to format properly (this was addressed in an edit to the question)

  • the mathematical formulas in the question are not clear in how they relate to the problem and are missing significant elements, e.g. alpha. This needs to be cleaned up. We have assumed that MPfonc gives the desired result and just simplified it.

The following corrects all the points and adds some minor improvements.

library(nls2)

set.seed(123) # for reproducibility

dens <- runif(100)
surf <- rnorm(100, 10, 2)
zone <- seq(1, 100, 1)
salaire <- rnorm(100, 1000, 3)

MPfonc <- function(alpha) {
  sapply(1:100, function(i) sum( (dens / (surf[i] - surf) ^ alpha)[-i], na.rm = TRUE ))
}

fo <- salaire ~ cbind(1, dens, surf, MPfonc(alpha))
gstart <- data.frame(alpha = c(-10, 10))

fm <- nls2(fo, start = gstart, alg = "plinear-random")

giving:

> fm
Nonlinear regression model
  model: salaire ~ cbind(1, dens, surf, MPfonc(alpha))
   data: parent.frame()
     alpha      .lin1  .lin.dens  .lin.surf      .lin4 
   0.90477 1001.20905   -0.50642   -0.12269    0.00681 
 residual sum-of-squares: 757.6

Number of iterations to convergence: 50 
Achieved convergence tolerance: NA

Note: Now that we have the starting value we can use it with nls like this:

nls(fo, start = coef(fm)["alpha"], alg = "plinear")

Update: Some code improvements, corrections and clarifications.

Upvotes: 3

Related Questions