Olá pessoal da lista! Estou lidando com um problema de certa forma parecido com o jogo SUDOKU. Resumindo nesse jogo tem uma área subdividida em quadriculas ou células que devem ser preenchidas seguindo condições de que os números de 1 a 9 não se repitam nas colunas e linhas. No meu caso, eu preciso interpolar espacialmente uma variável chamada de intensidade da chuva, entretanto, eu tenho a média dessa variável entre dois pontos que ligando eles formam uma linha, figura abaixo. A questão é a seguinte, como a partir da média eu obtenho os dados pontuais sobre as linhas que satisfaça a condição de que a média permaneça até mesmo nos pontos onde as linhas se cruzam. Ou seja, onde elas se cruzam precisa ser um valor que satisfaça a média da linha 1 e da linha 2. Depois disso preciso interpolar espacialmente esses pontos na área toda. Eu obtive uma solução, mas não sei se é a mais adequada e talvez para uma banco de dados maior ela pode ser impraticável. Resumindo, eu simulei valores nas linhas com a função grf() do geoR e depois com a função optim() variei as células (pixels) até que as condições sejam satisfeitas. Segue o código reproduzível e em anexo os arquivos referente ao shapefile que representa as linhas. proj="+proj=sterea +lat_0=52.15616055555555 +lon_0=5.38763888888889 +k=0.9999079 +x_0=155000 +y_0=463000 +ellps=bessel +units=km +no_defs " ## Colocar os arquivis anexados na mesma pasta da código R Links2l <- spTransform(readOGR("Links2l.shp"), CRS(proj)) Rmean <- c(2.9424, 4.3982) ## Criando uma imagem raster vazia ras <- raster(Links2l, res=0.1, ext=extent(Links2l)); ras[] <- rep(NA, ncell(ras)) ## Celulas (pixels) interceptadas pelas links cells_links <- cellFromLine(ras, Links2l) cells <- unique(unlist(cells_links)) ## Celulas compartilhadas entre os links cells_links[[1]][which(cells_links[[1]] %in% cells_links[[2]])] ## Simulando valores a partir da media nas celulas interceptadas pelas links set.seed(234) for(i in 1:nrow(Links_ams)){ ras[cells_links[[i]]] <- apply(grf(length(cells_links[[i]]), grid=xyFromCell(ras, cells_links[[i]]), nsim=10, mean = Rmean[i],cov.pars = c(3.7, 18.7), nugget = 0.37)$data, 1, mean) } plot(ras, col = gray.colors(10, start = 0.3, end = 0.9, gamma = 2.2, alpha = NULL)) lines(Links2l) ## Funcao de otimizacao opt <- function(par, Rmean){ res <- c() for(i in 1:length(Rmean)){ res[i] <- sum(mean(par[cells %in% cells_links[[i]]])-Rmean[i])^2 } sum(res) } ## Definir espaco de iteracao lower <- list() upper <- list() for(i in 1:nrow(Links_ams)){ lower[[i]] <- rep(mean(ras[cells_links[[i]]])-3.7^0.5, ncell(cells_links[[i]])) upper[[i]] <- rep(mean(ras[cells_links[[i]]])+3.7^0.5, ncell(cells_links[[i]])) } lower <- ifelse(unlist(lower)<0, 0, unlist(lower)) upper <- ifelse(unlist(upper)<0, mean(unlist(upper))+sd(unlist(upper)), unlist(upper)) ## Otimizar e atualizar as celulas para resolver a condicao cell_upd <- optim(par=ras[cells], fn=opt,lower = lower, upper=upper, Rmean=Rmean,method="L-BFGS-B")$par ras_old <- ras ras[cells] <- cell_upd ## Comparando celulas antigas com as novas, apenas a célula compartilhada entre os links nao é linear plot(ras_old[cells_links[[1]]], ras[cells_links[[1]]], ann=T) abline(a=0, b=1) plot(ras_old[cells_links[[2]]], ras[cells_links[[2]]], ann=T) abline(a=0, b=1) ## Verificando se as condicoes foram satisfeitas mean(ras[cells_links[[1]]]) Rmean[1] mean(ras[cells_links[[2]]]) Rmean[2] ## Resultado, valores iguais "SUDOKU" resolvido! ###################################################### ## Interpolacao das celulas pts <- rasterToPoints(ras) plot(pts) boxcox(pts~1) geo_link <- as.geodata(pts, coords.col = 1:2, data.col = 3) plot(geo_link, low =T) plot(geo_link, trend='1st', low =T) plot(geo_link, trend='2nd', low =T) vario <- variog(geo_link,max.dist=10,uvec=seq(0, 10, by=0.01)) plot(vario) link_fit <- list() link_fit$cte <- likfit(geo_link, cov.model="exp", trend ="cte",ini=c(0.5,5), nugget = 0.1) link_fit$lon <- likfit(geo_link, cov.model="exp", trend =geo_link$coords[, 1],ini=c(0.5,5), nug=0.1) link_fit$st <- likfit(geo_link, cov.model="exp", trend ="1st",ini=c(0.5,5), nug=0.1) link_fit$nd <- likfit(geo_link, cov.model="exp", trend ="2nd",ini=c(0.5,20), nug=0.1) sapply(link_fit, AIC) summary(link_fit$st) ## Predição na área krigagem grid <- xyFromCell(ras, 1:ncell(ras)) link_fit$st kr.link <- krige.conv(geo_link, loc=grid, krige=krige.control(obj=link_fit$st)) ras[] <- kr.link$predict image(ras) contour(ras, add=T) lines(Links2l, col=3, lwd=2) ## Krigagem ordinaria é BLUE, consequentemente condicoes sao satisfeitas mean(ras[cells_links[[1]]]) Rmean[1] mean(ras[cells_links[[2]]]) Rmean[2] Att. Abraco -- *Wagner Wolff, **PhD* "*Luiz de Queiroz**" College of Agriculture,* University of São Paulo Pádua Dias avenue11 | 13418-900| Piracicaba-SP| Brazil Phone: +55 19 982385582 <+55%2019%2098238-5582> http://orcid.org/0000-0003-3426-308X https://github.com/wwolff7 http://buscatextual.cnpq.br/buscatextual/visualizacv.do?id=K4463141A1
