¿Cómo realizar una regresión no lineal del pixel por pixel?

Question

El uso de nlsLM de Paquete 'minpack.lm' es sencillo para un ejemplo sencillo como esto:

MODEL:y=(exp(a*x+b*z+c)+d)^f
x=c(0.5,0.3,0.2,0.4)
z=c(0.1,0.6,1,0.9)
y=c(0.2,0.3,0.9,0.9)
fit=nlsLM(y~(exp(a*x+b*z+c)+d)^f,start = list(a = -.03, b = 0.5, c = 1,d=0.02,f=0.003))

Y no hay problema con eso. Mis datos reales son varios rásteres de un año.Así que tengo los datos de un año para x,y,z, pero como los rásteres que contienen píxeles,yo quiero hacer una regresión no lineal de un pixel por un pixel(regresión para cada uno de los píxeles no está relacionado con el siguiente píxel). datos de ejemplo:

   r <- raster(nrows=10, ncols=10); r <- setValues(r, 1:ncell(r))
   r1 <- raster(nrows=10, ncols=10);r1 <- setValues(r1, 1:ncell(r))
   r2 <- raster(nrows=10, ncols=10);r2 <- setValues(r1, 1:ncell(r))
st1=stack(r,r1,r2)

   rl <- raster(nrows=10, ncols=10); r <- setValues(r, 1:ncell(r))
   rS <- raster(nrows=10, ncols=10);r1 <- setValues(r1, 1:ncell(r))
   rT <- raster(nrows=10, ncols=10);r2 <- setValues(r1, 1:ncell(r))
st2=stack(rl,rS,rT)

  re <- raster(nrows=10, ncols=10); r <- setValues(r, 1:ncell(r))
  ru <- raster(nrows=10, ncols=10);r1 <- setValues(r1, 1:ncell(r))
  rg <- raster(nrows=10, ncols=10);r2 <- setValues(r1, 1:ncell(r))
st3=stack(re,ru,rg)

así que primero x sería el primer píxel de la st1(para todo el período temporal)

así que la primera y sería el primer píxel de st2(para todo el período temporal)

así que primero z sería el primer píxel de st3(para todo el período temporal)

hacer una regresión para este píxel, a continuación, hacer lo mismo para el resto de los píxeles

Como resultado puedo obtener un mapa (raster) 10*10 para cada uno de las mejores de parámetros(a,b,c,d,f)

Espero que quede claro y usted tiene cualquier ayuda!

Answer 1

Esto debe hacer lo que quiera.

El modelo que he utilizado no es el que se registra y probablemente no tiene sentido, sin embargo demostrar el principio.

# Set up the rasters
r1 <- r2 <- r3 <- r4 <- r5 <- r6 <- raster(nrows=10, ncols=10);
# Populate them with some values
r1 <- setValues(r1,runif(100,min=1,max=100));
r2 <- setValues(r2,runif(100,min=1,max=100));
r3 <- setValues(r3,runif(100,min=1,max=100));
r4 <- setValues(r4,runif(100,min=1,max=100));
r5 <- setValues(r5,runif(100,min=1,max=100));
r6 <- setValues(r6,runif(100,min=1,max=100));
# Stack them
st1 <- stack(r1,r2,r3);
st2 <- stack(r4,r5,r6);
# Set up the function
test <- function(r) {
    x <- r[1:3];
    y <- r[4:6];
    result <- c(NA,NA);
    try(result <- c(coef(nlsLM(y ~ a + b * x, start = list(a = 0.12345, b = 0.54321), na.action = na.omit))));
    result;
}
# Stack the stacks
s <- stack(st1,st2);
# Calculate a new Raster
rNLR <- calc(s, test);