Abstract
El uso de información georreferenciada se ha incrementado en los últimos años debido al incremento de la disponibilidad de este tipo de datos. Esta información georreferenciada o espacial tiene unas características particulares que requiere de una serie de herramientas estadísticas específicas. En este artículo haré una introducción a los distintos tipos de datos espaciales y cómo la estadística espacial nos puede ayudar a analizarlos. Además, comentaré algunas de las aplicaciones más importantes y los ilustraré con algunos ejemplos reales.
Keywords: Datos en retícula, Estadística espacial, Geoestadística, Procesos puntuales.
MSC Subject classifications: 62H11, 62M30.
Introducción
El geógrafo Waldo Tobler enunció en su primera ley de la geografía que “todas las cosas están relacionadas entre sí, pero las cosas más próximas están más relacionadas que las distantes”. La estadística espacial estudia procesos que ocurren en el espacio y de los que, de una manera u otra, se conoce su ubicación. Desde un punto de vista estadístico, esta dependencia basada en la distancia puede modelizarse de distintas maneras utilizando las herramientas proporcionadas por la estadística espacial. Tradicionalmente, suelen distinguirse tres áreas principales ((Cressie 2015)) dentro de la estadística espacial: la geoestadística, los procesos puntuales y los datos agregados.
La geoestadística se centra en el estudio de variables que cambian en el espacio de manera continua como, por ejemplo, la contaminación ambiental en una ciudad o la temperatura en la superficie del mar. En este caso, se dispone de observaciones en unos pocos puntos de medición, a partir de las cuales se intenta estimar la variable de interés en otros puntos de la región de estudio.
En el análisis de procesos puntuales se estudia la distribución de la ubicación de una serie de eventos, un patrón puntual, dentro de la región de estudio. Ejemplos típicos incluyen la distribución de especies de árboles o la de casos de una enfermedad en una ciudad.
Por último, el análisis de datos agregados se refiere al estudio de variables que no están necesariamente relacionadas con una ubicación en el espacio sino que están asociadas a determinadas divisiones del mismo. Por ejemplo, estaríamos hablando del estudio del número de semillas que germinan en cada una de las parcelas en las que se ha dividido un campo más extenso o del número de casos de una enfermedad por barrios u otra división administrativa similar.
Geoestadística
La geoestadística se centra en el estudio de variables espaciales que cambian de manera continua dentro de la región de estudio. Para abordar este tipo de problemas se tienen, normalmente, observaciones de las variables de interés en unas pocas localizaciones del espacio. A partir de estas observaciones es posible estimar el valor de la variable de interés en cualquier punto de la región de estudio.
Los primeros métodos geoestadísticos se desarrollaron dentro del contexto de la minería en la década de 1960, por el ingeniero sudafricano D. G. Krige, con el objetivo de estimar las reservas mineras. Estas ideas fueron formalizadas por G. Matheron en lo que se conoce comúnmente como kriging ((Chilès and Desassis 2018)).
Si consideramos un proceso estocástico \(Z(x)\), donde \(x\) indica la ubicación del punto del espacio \(D\) en el que se observó la variable, es común asumir que el proceso es estacionario de segundo orden, es decir, que en toda la región de estudio \(D\) tiene media \(\mu\) y varianza \(\sigma^2\) constantes. Además, suele considerarse también que el proceso es isotrópico, lo que implica que la covarianza entre dos observaciones cualesquiera solo depende de la distancia que las separa (y no de la dirección o su posición relativa).
El estimador del kriging ordinario se basa en una suma ponderada de los valores observados en \(n\) puntos: \[\hat{Z}(x) = \sum_{i=1}^n \lambda_i Z(x_i).\qquad(1)\]Los pesos \(\{\lambda_i\}_{i=1}^n\) son tales que su suma es uno y se obtienen de manera que se minimice la varianza del error de estimación. Estos pesos dependen de lo que se conoce como (semi)variograma, que mide la dependencia entre dos puntos \(x_1\) y \(x_2\) cualesquiera del espacio: \[\gamma(x_1, x_2) = \frac{1}{2}E[(Z(x_1) – Z(x_2))^2]\qquad(2)\]donde \(E[\cdot]\) denota la esperanza. En un proceso estacionario e isotrópico el variograma solo depende de la distancia entre dos puntos.
Las ideas detrás del kriging ordinario pueden generalizarse para considerar procesos en los que la media no es constante y también para considerar procesos anisotrópicos, en los que la posición relativa entre dos puntos, y no solo su distancia, influye en la dependencia espacial de los mismos. Por ejemplo, este sería el caso de industrias contaminantes en las que el viento esparce la contaminación en una determinada dirección.
A la hora de construir modelos geoestadísticos es muchas veces indispensable contruir un modelo multivariante para estudiar tanto la tendencia (media) como la depdencia (covarianza) entre las distintas observaciones en función de la distancia entre dos puntos. Los procesos gaussianos son modelos que nos van a permitir analizar datos geoestadísticos de una manera muy flexible.
Diremos que un proceso estocástico definido sobre una región \(D\) es un proceso gaussiano cuando, si tenemos \(n\) observaciones del proceso en otras tantas localizaciones \(\{x_i\}_{i=1}^n\), el vector de observaciones \((y_1,\ldots,y_n)\) sigue una distribución multivariante de media \(\mu\) y varianza \(\Sigma\). Tanto la media como la varianza pueden modelizarse de distintas maneras para incluir dependencia espacial. La variación a gran escala o tendencia suele incluirse como parte de la media mientras que el término del error espacial suele modelizarse utilizando la varianza \(\Sigma\) ((Cressie 2015)).
La media del proceso gaussiano se puede modelizar utilizando la ubicación de las observaciones, así como variables topográficas, ambientales o de otro tipo. Para modelizar la varianza del proceso gaussiano es común utilizar una función de covarianza. Para procesos estacionarios e isotrópicos, la función de covarianza de Matérn ((Krainski et al. 2018)) ofrece un marco bastante conveniente dada su flexibilidad a la hora de modelizar distintos tipos de dependencia espacial.
De manera más general, estos modelos pueden extenderse para considerar distribuciones no gaussianas para los datos observados. Para ello se puede recurrir a los modelos lineales generalizados ((Diggle and Ribeiro 2007)) que permiten asumir distintas distribuciones para la variable respuesta y que modelizan la media de esta distribución en función de un predictor lineal \(\eta_i\) que puede incluir términos geoestadísticos: \[\eta_i = X_i^{\top} \beta + S(x_i) + \ldots\qquad(3)\]En la ecuación anterior \(X_i\) es un vector de covariables medidas en el punto \(x_i\) y \(S(\cdot)\) es un efecto espacial que varía continuamente en el espacio y que puede ser modelizado utilizando términos estadísticos.
Estudio de la contaminación ambiental
Un ejemplo clásico del uso de modelos geoestadísticos es la estimación de la contaminación a partir de una serie de mediciones en unos pocos puntos. Un ejemplo es el estudio de la contaminación por metales pesados cerca del pueblo de Stein (Limburgo, Países Bajos) que está disponible en el paquete gstat ((Pebesma 2004)) de \(\texttt{R}\) ((R Core Team 2021)). La Figura 1 muestra la ubicación de los puntos de muestro, que se encuentran en las inmediaciones del río Mosa. En estos puntos se tomaron muestras del suelo y se midió la concentración de varios metales pesados (cadmio, cobre, zinc y plomo). Además, se obtuvieron otras variables como la altitud del terreno, la distancia al río y la frecuencia con la que se inundan los terrenos.
En la Figura 1 se puede apreciar un mayor nivel de concentración de zinc (en ppm) en los puntos más cercanos al río, que también son los que se inundan con mayor frecuencia.
Mediante las técnicas de modelización geoestadística es posible obtener una estimación de la concentración de zinc en una malla regular sobre la región de estudio. En primer lugar, el variograma nos permitirá ver cómo es la dependencia espacial de la variable de interés que, en este caso, es el logartimo de la concentración de zinc (en ppm). La Figura 2 muestra el variograma. Los puntos muestran el variograma empírico (obtenido a partir de los datos) mientras que la línea negra muestra un ajuste a un variograma esférico. A partir de los 900 metros de distancia el variograma se estabiliza, lo que indica que a partir de esa distancia no existe dependencia espacial entre los puntos.
La Figura 3 muestra los valores estimados usando kriging ordinario junto con el error típico de las estimaciones (que suele incrementarse al alejarnos de los puntos de muestreo).
Epidemiología espacial
El uso de modelos geoestadísticos en epidemiología es muy interesante porque permite modelizar que el riesgo de sufrir una determinada enfermedad varía de manera continua en la región de estudio. Por ejemplo, para el estudio de la variación espacial de la incidencia de enfermedades tropicales desatendidas ((Diggle and Ribeiro 2007)) se utilizan modelos en los que la variable respuesta es un conteo de casos que depende de varios efectos, entre ellos un patrón espacial estimado con modelos geoestadísticos.
La Figura 4 muestra la prevalencia de la malaria (es decir, la proporción de niños con malaria en sangre) en las distintas aldeas en las que se tomaron los datos, que se obtuvieron durante un estudio sobre la prevalencia de malaria en niños en aldeas de Gambia ((Diggle et al. 2002)). Ahí se puede ver un patrón espacial bastante claro, con una zona de alta prevalencia al este del país.
En este caso se modelizó la variable \(Y_{ij}\), que representa si un niño \(j\) de una aldea \(i\) tiene malaria (\(Y_{ij}=1\)) o no (\(Y_{ij} =0\)), utilizando una distribución de Bernoulli con probabilidad \(p_{ij}\): \[\begin{array}{rcl}
Y_{ij} & \sim & Be(p_{ij})\\
\log(p_{ij}/ (1 – p_{ij})) & = & X_{ij}^{\top} \beta + U_i + S(x_i)
\end{array}\qquad(4)\]donde \(X_{ij}\) es un vector de covariables con coeficientes, \(\beta\) \(U_i\) representa un término de errores independientes y \(S(x_i)\) un efecto espcial que varía de manera continua. De esta manera, se tienen en cuenta los factores de riesgo individuales de cada niño (mediante las covariables \(X_{ij}\), la probable correlación espacial en los datos (mediante el efecto espacial \(S(x_i)\)) y los factores intrínsecos a las aldeas (con el efecto \(U_i\)).
Modelos con barreras
Los modelos geoestadísticos más sencillos consideran que la variación de la variable de interés solo va a depender de la distancia entre dos puntos y, posiblemente, de su posición relativa. Sin embargo, determinados procesos espaciales pueden estar afectados por barreras dentro de la región de estudio. Este es el caso cuando se estudia el nivel de ruido en una ciudad, por ejemplo. Un modelo geoestadístico usual solamente tendría en cuenta la distancia a los puntos de muestreo a la hora de estimar el nivel de ruido, sin ignorar que existen barreras arquitectónicas que afectan a la propagación del ruido (como, por ejemplo, los edificios). (Bakka et al. 2019) han propuesto una metodología para modelizar datos geoestadísticos y tener en cuenta las barreras de distinto tipo. Estos modelos han sido utilizados para el análisis de la distribución espacial de cetáceos ((Martínez-Minaya et al. 2019)) donde las islas de un archipiélago presentan barreras. También han sido utilizados para estimar el nivel de ruido en el centro de las ciudades ((Krainski et al. 2018)), donde los edificios actúan como barreras.
La Figura 5 muestra la correlación espacial de un modelo basado en un proceso espacial basado en modelos con barreas y un modelo estacionario calculada sobre las calles de la zona centro de la ciudad de Albacete. Como puede verse, esta última no tiene en cuenta la presencia de obstáculos, como los edificios (en gris), a la hora de obtener la dependencia espacial entre dos puntos.
Procesos puntuales
El análisis de procesos puntuales ((Diggle 2013)) se centra en el estudio de la distribución espacial de la ocurrencia de determinados eventos. Ejemplos típicos incluyen la ubicación de árboles en un bosque o la de casos de cáncer en una ciudad.
En general, las preguntas que surgen al analizar un proceso puntual se centran en si la aparición de los sucesos es uniforme en el espacio o si existen zonas en las que los sucesos tienden a ocurrir con mayor o menor frecuencia. En el caso más sencillo, se asume que los eventos ocurren independientemente unos de otros y que su distribución en toda la región de estudio \(D\) es uniforme, de manera que el número de eventos observado sigue una distribution de Poisson de media \(\lambda\). Al parámetro \(\lambda\) se le llama la intensidad del proceso puntual, y representa el número de casos por unidad de área. A un proceso que cumple estas propiedades se le conoce como proceso homogéneo de Poisson con intensidad \(\lambda\).
La extensión natural al proceso homogéneo de Poisson es permitir que la intensidad no sea constante y que haya zonas en las que se observen más eventos que otras, de manera que la intensidad dependa ahora de la ubicación de los puntos, es decir, ahora la intensidad es \(\lambda(x)\) con \(x\) un punto de la región de estudio \(D\). En este proceso se pueden observar agregaciones de puntos, pero esta agregación se debe únicamente al caracter no homogéneo de la intensidad y en ningún casos a la existencia de dependencia entre los eventos, que vuelven a aparecer de manera independiente.
La Figura 6 muestra un ejemplo de un proceso puntual de la distribución espacial de tres especies de árboles en Lansing Woods, Clinton County, Michigan (EE.UU.). Estos datos se han obtenido del banco de datos lansing del paquete spatstat ((Baddeley, Rubak, and Turner 2015)) del software estadístico R ((R Core Team 2021)). Como puede observarse, la distribución de los árboles no parece ser homogénea e, incluso, parece existir un patrón espacial distinto para el roble negro y el arce.
Además de las marcas que indican el tipo de evento al que hacen referencia (por ejemplo, distintas especies de árboles) ((Gomez-Rubio, Cameletti, and Finazzi 2015)), un proceso puntual puede contener también información adicional sobre los eventos (por ejemplo, altura del arbol). Esta información se puede incorporar a la modelización de la intensidad de manera que obtengamos \[\log(\lambda(x)) = S_0(x) + \beta_0 + \beta_1 z(x) + \ldots\qquad(5)\] En la expresión anterior, \(S_0(x)\) representa un efecto espacial que varía de forma continua, \(\beta_0\) la intersección del modelo y \(\beta_1\) es el coeficiente de una covariable \(z(x)\) que puede variar espacialmente. La estimación de estos efectos puede hacerse de varias maneras.
Por ejemplo, el efecto espacial \(S_0(x)\) representa una superficie suave que puede estimarse de manera paramétrica o no paramétrica mediante, por ejemplo, una suavización kernel ((Diggle 1985)). La Figura 7 muestra la estimación de la intensidad utilizando estimación kernel. En esta gráfica puede apreciarse mejor la distribución no uniforme de las tres especies, además de las diferencias entre el arce y las otras dos especies de árboles.
Ecología
Un problema interesante en ecología es estudiar la distribución de especies animales y vegetales. En el caso de estas últimas, existen registros extensos (como el Inventario Forestal Nacional) con la ubicación exacta de árboles de distintas especies vegetales. De esta manera, se puede estudiar la distribución por separado de estas especies vegetales y, además, si existe algún tipo de interacción entre ellos o no.
Es posible que distintas especies vegetales compitan por el espacio, lo que origine un efecto de repulsión entre ellos, de manera que sea poco frecuente que aparezcan en las mismas zonas. También puede darse el caso de que distintas especies tengan una relación de tipo simbiótico que haga que estas especies aparezcan de manera conjunta frecuentemente.
La información sobre la ubicación de especies vegetales se suele completar con otro tipo de información sobre el terreno obtenida de distintas maneras. Por ejemplo, no es difícil incluir información sobre el terreno (altitud, pendiente, etc.) que permita enriquecer el banco de datos. De esta manera, es posible estimar el efecto del terreno en la distribución espacial de especies vegetales.
La función \(K\) de Ripley ((Baddeley, Møller, and Waagepetersen 2000) y (Ripley 1977)) permite detectar si un patrón está más agregado (o disperso) con respecto a un modelo de Poisson con intensidad \(\lambda(x)\). Esto permite contrastar si existen mecanismos de interacción entre los eventos, lo que requeriría emplear modelos más complejos para su análisis ((Baddeley, Rubak, and Turner 2015), (Diggle 2013) y (Illian et al. 2008)).
Epidemiología
Una de las aplicaciones más interesantes de los procesos puntuales aparece en epidemiología. Cada vez es más común tener información detallada e individualizada sobre enfermos o fallecidos por una determinada causa. De esta manera, es posible conocer el lugar de residencia de enfermos y posibles factores de riesgo. Esta información se puede complementar con información socioeconómica proporcionada por el Instituto Nacional de Estadística por distritos censales. Además, España mantiene un registro de información de industrias contaminantes, que incluye la ubicación, tipo de emisión y cantidad emitida, y que se puede consultar en https://prtr-es.es.
Así, es posible estudiar la distribución espacial de los casos de, por ejemplo, cáncer teniendo en cuenta distintos factores de riesgo socioeconómicos y la distancia a posible industrias contaminantes (que se usa frecuentemente como una aproximación a la exposición).
Medicina
Como último ejemplo, comentaremos el uso de procesos puntuales en medicina. Se trata, sobre todo, de aplicaciones en imagen médica, a partir de la cual se determina la úbicación de determinadas célculas en un tipo de tejido. Estudiar la distribución de estas células puede ayudar a determinar mecanismos fisiológicos subyacentes. Diversos autores ((Baddeley, Rubak, and Turner 2015), (Diggle 2013) y (Illian et al. 2008)) proporcionan numerosos ejemplos de estas aplicaciones.
Datos agregados
En el análisis de procesos puntuales el proceso espacial se estudia de manera individualizada ya que se conoce la ubicación de todos los eventos de interés. Otras veces es imposible acceder a esta información desagregada y los datos espaciales aparecen agregados de acuerdo a una división de la región de estudio. Un caso típico es considerar un terreno dividido en varias parcelas, en las que se cuenta el número de árboles de una determinada especie.
El Instituto Nacional de Estadística proporciona información espacial agregada a distintos niveles administrativos dependiendo del tipo de datos. Por ejemplo, estadísticos resumen de muchas variables socioeconómicas obtenidas en el censo están disponibles por secciones censales. También existen datos de población y defunciones (desagregadas por edad, sexo y causa de defunción) a nivel municipal.
En general, el análisis de datos agregados se basa en distintos de modelos de regresión en los que, de distintas maneras, se modeliza la correlación en la variable respuesta entre regiones vecinas. Este concepto de vecindad entre regiones se puede definir de varias formas, aunque es muy común considerar que dos regiones son vecinas si sus límites se tocan en, al menos, uno o dos puntos.
Esta relación de vecindad se puede representar mediante un grafo, tal y como se muestra en la Figura 8. En ella se han representado los barrios de la ciudad de Columbus (Ohio, EE.UU.) y se han unido con arcos los centroides de los barrios que son vecinos. Estos datos se han obtenido del paquete spData ((Bivand, Nowosad, and Lovelace 2019)) de R.
Los modelos que se utilizan para el análisis de datos agregados puede derivarse del modelo lineal general en el que se incluye un término de error de manera que los errores están correlacionados entre área vecinas. Cuando la variable respuesta no sigue una distribución normal, es posible utilizar modelos lineales generalizados ((McCullagh and Nelder 2019)), por lo que es posible considerar variables respuesta binarias o que representen un conteo de eventos. De esta manera, es posible estudiar el efecto de las covariables incluidas en los modelos y también el efecto espacial al incluir esa correlación entre vecinos.
Un modelo lineal mixto que incluya una componente espacial podría estar definido mediante una distribución normal multivariante de media \(\mu\) y matriz de varianzas \(\sigma^2 I\) (con \(\sigma^2\) un parámetro de escala e \(I\) la matrix identidad). De esta manera, la media podría incluir efectos espaciales así: \[\mu = X \beta + u\qquad(6)\]con \(X\) una matriz de covariables, \(\beta\) sus coeficientes y \(u\) un vector de efectos espaciales. Este vector \(u\) suele estar definido como una normal multivariante de media cero y varianza \(\Sigma\). La dependencia espacial se incluye a la hora de definir la matriz \(\Sigma\), que suele depender de una matriz de pesos espaciales \(W\). Por ejemplo, tomando \(\Sigma\) como \((I – \rho W)^{-1}\) se incluye cierta dependencia espacial que depende del parámetro \(\rho\). En el caso de que \(\rho\) sea cero tendríamos un efecto sin ningún tipo de dependencia espacial.
Econometría espacial
Una de las áreas en las que más se han utilizado los modelos espaciales de datos agreagados es la econometría. El desarrollo de modelos de econometría espacial ((Anselin 1988) y (LeSage and Pace 2009)) ha ido en paralelo al desarrollo de modelos espaciales dentro del mundo de la estadística, tal vez por la necesidad de modelizar explícitamente dentro del modelo la dependencia espacial con regiones vecinas.
Los modelos de econometría espacial pueden incluir un término autorregreseivo sobre el término del error o la variable respuesta. En el primer caso, el modelo se escribiría (en forma matricial) así: \[Y = X \beta + u; \ u = \lambda W u + e\qquad(7)\]donde \(Y\) representa la variable respuesta, \(u\) es el efecto autocorrelado y \(e\) es un efecto aleatorio de media cero y varianza igual a \(\sigma^2_e I\) (es decir, errores independientes con varianza \(\sigma^2_e I\)). Este modelo se puede reescribir tal que \(u\) es un error normal multivariante de media cero y varianza \(\sigma^2_e I((I – \lambda W^{\top})(I – \lambda W))^{-1}\).
Cuando el término autorregresivo es en la variable respuesta, el modelo es \[Y = \rho W Y + X \beta + e .\qquad(8)\]De esta manera, el modelo se puede reescribir así: \[Y = (I – \rho W)^{-1}X \beta + u\qquad(9)\]donde \(u\) es un error autorregresivo como el definido arriba.
Este modelo presenta la característica de que el término de las covariables va multiplicado por \((I – \rho W)^{-1}\), de manera que el efecto de las covariables no solo afectaría al área en la que se han medido, sino también a los vecinos.
Esta es una característica que tienen los modelos clásicos de econometría espacial, puesto que permiten estimar el efecto que tendría el cambio en una covariable en la variable respuesta en las regiones vecinas. Este efecto se conoce como efecto derrame (spillover effect, en inglés).
Atlas de mortalidad
Una aplicación muy interesante del análisis espacial de datos agregados es el desarrollo de atlas de mortalidad dentro de la epidemiología espacial ((Martínez-Beneito and Botella-Rocamora 2019)). A partir de la información de población y mortalidad, es posible estimar el riesgo relativo en cada división administrativa utilizando estos modelos espaciales. Es común que estos modelos consideren dos términos de error: uno con errores correlacionados espacialmente y otro con errores independientes entre regiones ((Besag, York, and Mollié 1991)) (como ya vimos en el caso de modelos geoestadísticos en epidemiología). De esta manera, es posible tener en cuenta que regiones vecinas tengan riesgos similares, y a la vez permitir diferencias intrínsecas entre regiones.
En el desarrollo de atlas de mortalidad se suelen estudiar distintas causas de muerte. Esto permite comparar los patrones espaciales de distintas causas y determinar si existen patrones similares. Enfermedades con factores de riesgo similares suelen mostrar patrones espaciales similares.
En concreto, para cada área podemos tener el número de casos observado \(O_i\), el de casos esperados \(E_i\) y, posiblemente, un vector de covariables \(X_i\). El modelo, al ser un conteo, suele considerar que los casos observados siguen una distribución de Poisson con media \(E_i\theta_i\), con \(\theta_i\) el riesgo relativo. Este riesgo relativo se puede modelizar en función de las covariables y otros efectos así: \[\log(\theta_i) = X_i \beta + u_i + v_i\qquad(10)\]donde \(u_i\) es un efecto espacial y \(v_i\) efectos independientes. De esta manera, \(u_i\) modeliza la variación espacial del riesgo mientras que \(v_i\) tiene en cuenta las características intrínsecas de cada area. Por otro lado, el vector de coeficientes \(\beta\) permite medir el efecto de los factores de riesgo incluidos en el vector de covariables \(X_i\).
La Figura 9 muestra un ejemplo del atlas de mortalidad desarrollado por investigadores de la Fundación para el Fomento de la Investigación Biomédica de la Comunidad Valenciana (FISABIO) en el que se muestra el riesgo relativo suavizado considerando todas las causes de muerte en España a nivel municipal en el período 1989-2014. El atlas completo incluye análisis de varias causas de muerte por sexos y puede consultarse en https://medea3.shinyapps.io/atlas_nacional/.
Estos modelos se pueden extender para considerar efectos temporales y también es posible desarrollar modelos espaciales multivariantes ((Botella-Rocamora, Martinez-Beneito, and Banerjee 2015)) en los que además de la distribución espacial del riesgo se modeliza una posible correlación entre distintas causas de muerte, lo que permite cuantificar su posible relación.
Conclusiones
Como hemos visto, el análisis de datos espaciales o georreferenciados puede hacerse de distinta manera en función del problema y tipos de datos disponibles. La estadística espacial proporciona una amplia cantidad de métodos de análisis para tener en cuenta la probable dependencia espacial entre las observaciones. Además del desarrollo metodológico, herramientas como el lenguaje de programación R ((R Core Team 2021)) han permitido realizar análisis de datos espaciales cada vez más complejos.
Agradecimientos
Este trabajo ha sido realizado gracias al apoyo del proyecto PID2019-106341GB-I00 financiado por el Ministerio de Ciencia e Innovación (España) y una ayuda a grupos de investigación de la Universidad de Castilla-La Mancha.