Una guida passo-passo dei principali comandi per l'analisi statistica dei dati spaziali

spatial_analysis.R

# ===============================================
# === ANALISI SPAZIALE CON R                  ===
# === AUTHORS: Edbali Ossama, Gazzini Thommy  ===
# ===============================================

# =========================================
# === PARTE 1: VISUALIZZAZIONE DATASETS ===
# =========================================

# Carica dataset "meuse" (fiume Mosa)
data(meuse)

# Scatterplot di "meuse" (non leggibile)
plot(meuse)

# Converto il dataframe nella classe SpatialPointsDataFrame --> visualizzazione migliore
coordinates(meuse) <- c('x', 'y')
# Controllo classe attuale del dataset
class(meuse)
# Visualizzazione di informazioni sulle coordinate
plot(meuse)
# Aggiunta di altri attributi con spplot
names(meuse)
# Distribuzione di zinco
spplot(meuse, 'zinc')
# Grafico con scala logaritmica (log scale)
spplot(meuse, 'zinc', do.log = TRUE)

# Visualizzazione alternativa: 'bubble plot'
bubble(meuse, 'zinc')

# Ora definiamo "il contorno" del fiume
data(meuse.riv)
meuse.lst <- list(Polygons(list(Polygon(meuse.riv)), 'meuse.riv'))
meuse.sr <- SpatialPolygons(meuse.lst)
plot(meuse.sr, col='grey')

# Usando una 'grid' (grid plot)
data(meuse.grid)
coordinates(meuse.grid) <- c('x', 'y')
meuse.grid <- as(meuse.grid, "SpatialPixels")
image(meuse.grid, col='grey')
title('grid')

# Ora mettiamo tutte i tipi di visualizzazione assieme
image(meuse.grid, col='grey')
plot(meuse.sr, col = 'grey', add=TRUE)
plot(meuse, add=TRUE)

# ========================
# === ANALISI SPAZIALE ===
# ========================

# Uso del pacchetto spatstat per l'analisi spaziale di 
# pattern di punti, model-fitting e test
data(swedishpines)
X <- swedishpines
plot(X)
# Informazioni sul dataset
X
summary(X)
# Info sull'intensità dei punti
# Il valore 10 rappresenta la deviazione
# standard usata dalla matrice che
# effettua il filtro gaussiano
plot(density(X, 10))

# Quadrat counting
Q <- quadratcount(X, nx=4, ny=3)
plot(X)
plot(Q, add=TRUE, cex=2)

# K function
plot(Kest(X))

# Quadrat counting: seconda version
data(bei)
Z <- bei.extra$grad
b <- quantile(Z, probs=(0:4)/4)
Zcut <- cut(Z, breaks=b, labels=1:4)
V <- tess(image=Zcut)
plot(V)
plot(bei, add=TRUE, pch="+")

qb <- quadratcount(bei, tess=V)
plot(qb)

# Simulazione di un processo di Poisson
plot(rpoispp(100))

# Test del chi quadro di CSR usando i quadrati di campionamento
data(nztrees)
M <- quadrat.test(nztrees, nx=3, ny=2)
plot(nztrees)
plot(M, add=TRUE, cex=2)

# Formal check
data(bei)
fit <- ppm(bei,trend=~x)
M <- quadrat.test(fit, nx=4, ny=2)
M

# Informal check
plot(bei, pch=".")
plot(M, add=TRUE, cex=1.5, col="red")

# Metodo di Monte Carlo e inviluppi della funzione K
data(cells)
E <- envelope(cells, Kest, nsim=39, rank=1)
E
plot(E, main="Inviluppi")