Pessoal, acho que não formulaei bem o problema. Por favor, considerem a mensagem abaixo como a correta. Agradeço se alguém puder dar alguma indicação. Abraços. library(raster) library(logistf) # Criando 2 rasters para reproduzir o que estou enfrentando r <- raster(nrow=10, ncol=10) # Variável-raster resposta (binária) com 9 bandas/layers s1 <- lapply(1:9, function(i) setValues(r, sample(0:1,ncell(r),replace = T))) s1 <- stack(s1) # um raster-variável preditora val <- sample(0:60,ncell(r),replace = T) s2 <- raster(nrow=10, ncol=10,vals=val) plot(s1) plot(s2) # Uma segunda variável preditora armazenada como vetor numérico: 9 valores assim como existem 9 bandas/layres em s1 (não é coincidência) exp_2 <- c(27.00,30.02,31.07,32.72,33.73,35.12,35.65,36.06,38.32) Agora, eu quero ajustar o modelo usando o logistf (eu tenho motivos estatísticos para isso, não reproduzido aqui) usando a função calc do pacote raster. É aí que o mistério mora. O raciocínio na formulação do problema é que cada célula/pixel seguiria a segunite fórmula glm: s1 ~ 27.00 + pixel correspondente em s2 + 27.00:pixel correspondente em s2 s1 ~ 30.02 + pixel correspondente em s2 + 30.02:pixel correspondente em s2 s1 ~ 31.07 + pixel correspondente em s2 + 31.07:pixel correspondente em s2 ... e assim por diante para todos os 9 valores de exp_2. Eu tentei algo assim: fun <- function (x) { logistf (x ~ exp_2 + s2 + exp_2:s2)$coeficients } coefs <- calc(s1, fun) Mas estava claro que não funcionaria. Provavelmente por causa das diferentes dimensões/número de rasterlayers neste conjunto de dados. Estou recebendo a mensagem de erro: Error in .calcTest(x[1:5], fun, na.rm, forcefun, forceapply) : cannot use this function *Jefferson Ferreira-Ferreira, **PhD (abd)* *Geographer* *Ecosystem Dynamics Observatory <http://tscanada.wix.com/ecodyn> - EcoDyn/UNESP* *Department of **Geography * *Institute of Geosciences and Exact Sciences** (IGCE)* *São Paulo State University (UNESP)* *Rio Claro, SP - Brazil* Em 9 de agosto de 2017 16:07, Jefferson Ferreira-Ferreira <jecogeo@gmail.com

escreveu:

Pessoal,

Estou tentando ajustar a Estimativa de Máxima Máxima Penalizada de Firth para a regressão logística implementada no pacote logistf para um conjunto de rasters (imagens), mas estou com problemas para implementá-lo usando rasters e vetores numéricos. Para ilustrar o problema:

library(raster) library(logistf)

# Criando 2 rasters para reproduzir o que estou enfrentando r <- raster(nrow=10, ncol=10)

# Variável-raster resposta (binária) com 9 bandas/layers s1 <- lapply(1:9, function(i) setValues(r, sample(0:1,ncell(r),replace = T))) s1 <- stack(s1)

# um raster-variável preditora val <- sample(0:60,ncell(r),replace = T) s2 <- raster(nrow=10, ncol=10,vals=val)

plot(s1) plot(s2)

# Uma segunda variável preditora armazenada como vetor numérico: 9 valores assim como existem 9 bandas/layres em s1 (não é coincidência) exp_2 <- c(27.00,30.02,31.07,32.72,33.73,35.12,35.65,36.06,38.32)

Agora, eu quero ajustar o modelo usando o logistf (eu tenho motivos estatísticos para isso, não reproduzido aqui) usando a função calc do pacote raster. É aí que o mistério mora.

O raciocínio na formulação do problema é que cada célula/pixel seguiria a segunite fórmula glm:

s1/layer1 ~ 27.00 + pixel correspondente em s2 + 27.00:pixel correspondente em s2 s1/layer2 ~ 30.02 + pixel correspondente em s2 + 30.02:pixel correspondente em s2 s1/layer3 ~ 31.07 + pixel correspondente em s2 + 31.07:pixel correspondente em s2

... e assim por diante para todas as 9 bandas/layers do meu raster-variável resposta s1, que são pareadas com valores de exp_2.

Eu tentei algo assim:

fun <- function (x) { logistf (x ~ exp_2 + s2 + exp_2:s2)$coeficients } coefs <- calc(s1, fun)

Mas estava claro que não funcionaria. A parte complicada é infomar ao script que eu quero que cada valor de exp_2 seja usado com cada camada raster de s1 para este modelo.

Qualquer idéia será muito útil. Idéias?

Obrigado e abraços.

*Jefferson Ferreira-Ferreira, **PhD (abd)*

*Geographer*

*Ecosystem Dynamics Observatory <http://tscanada.wix.com/ecodyn> - EcoDyn/UNESP* *Department of **Geography * *Institute of Geosciences and Exact Sciences** (IGCE)* *São Paulo State University (UNESP)* *Rio Claro, SP - Brazil*