Pada tugas ini, diminta untuk:
Melakukan analisis terhadap data titik (boleh dari R package atau sumber lain)
Menjelaskan hasil analisis (minimal mencakup metode kuadran dan fungsi-K)
Menuliskannya dan mengumpulkan dalam bentuk html, lengkap dengan keterangannya.
Dataset nztrees yang tersedia di paket
spatstat merupakan data pola titik (planar point
pattern) yang merekam lokasi 64 pohon black beech
(Nothofagus solandri) di sebuah petak hutan alami di New
Zealand. Data ini disajikan dalam bentuk koordinat dua dimensi (x, y)
dengan satuan meter, yang terletak pada area berbatas persegi panjang
berukuran sekitar 56,5 × 38,5 meter. Secara umum, dataset ini banyak
digunakan sebagai contoh klasik dalam analisis spasial untuk mengkaji
pola distribusi pohon, apakah tersebar secara acak, cenderung
mengelompok, atau justru teratur. Melalui pendekatan analisis seperti
metode kuadran maupun fungsi-K, pola sebaran pohon dalam dataset ini
dapat dievaluasi sehingga memberi gambaran tentang dinamika ekologis dan
interaksi spasial antarindividu pohon di dalam hutan tersebut.
## Warning: package 'spatstat' was built under R version 4.5.1
## Loading required package: spatstat.data
## Warning: package 'spatstat.data' was built under R version 4.5.1
## Loading required package: spatstat.univar
## Warning: package 'spatstat.univar' was built under R version 4.5.1
## spatstat.univar 3.1-4
## Loading required package: spatstat.geom
## Warning: package 'spatstat.geom' was built under R version 4.5.1
## spatstat.geom 3.5-0
## Loading required package: spatstat.random
## Warning: package 'spatstat.random' was built under R version 4.5.1
## spatstat.random 3.4-1
## Loading required package: spatstat.explore
## Warning: package 'spatstat.explore' was built under R version 4.5.1
## Loading required package: nlme
## spatstat.explore 3.5-2
## Loading required package: spatstat.model
## Warning: package 'spatstat.model' was built under R version 4.5.1
## Loading required package: rpart
## spatstat.model 3.4-0
## Loading required package: spatstat.linnet
## Warning: package 'spatstat.linnet' was built under R version 4.5.1
## spatstat.linnet 3.3-1
##
## spatstat 3.4-0
## For an introduction to spatstat, type 'beginner'
## Planar point pattern: 86 points
## window: rectangle = [0, 153] x [0, 95] feet
## Planar point pattern: 86 points
## Average intensity 0.005916753 points per square foot
##
## Coordinates are integers
## i.e. rounded to the nearest foot
##
## Window: rectangle = [0, 153] x [0, 95] feet
## Window area = 14535 square feet
## Unit of length: 1 foot
Metode kuadran digunakan untuk menganalisis pola sebaran titik dengan cara membagi area studi ke dalam beberapa sel (kuadran) yang ukurannya sama. Langkah-langkahnya sebagai berikut:
Interpretasi VMR
Hipotesis
Uji Statistik
Jika \(m < 30\), gunakan pendekatan Chi-Square: \[ \chi^2 = \frac{(m-1) S^2}{\bar{x}} = (m-1) \times VMR \]
Jika \(m \geq 30\), gunakan pendekatan Z: \[ Z = \sqrt{\frac{m-1}{2}} \; (VMR - 1) \]
Keputusan (taraf nyata \(\alpha = 5\%\))
## x
## y [0,38.2] (38.2,76.5] (76.5,115] (115,153]
## (63.3,95] 4 5 4 12
## (31.7,63.3] 8 8 8 10
## [0,31.7] 8 7 5 7
## [1] 0.807611
##
## Chi-squared test of CSR using quadrat counts
##
## data:
## X2 = 8.8837, df = 11, p-value = 0.7348
## alternative hypothesis: two.sided
##
## Quadrats: 4 by 3 grid of tiles
##
## Chi-squared test of CSR using quadrat counts
##
## data:
## X2 = 8.8837, df = 11, p-value = 0.3674
## alternative hypothesis: regular
##
## Quadrats: 4 by 3 grid of tiles
Berdasarkan metode kuadran dengan pendekatan uji Khi-Kuadrat, sebaran
titik pada data nztrees tidak berbeda signifikan
dari pola acak (CSR). Dengan kata lain, distribusi pohon di
area studi ini dapat dianggap acak.
Secara matematis, fungsi \(K\) empiris didefinisikan sebagai:
\[ \hat{K}(r) = \frac{|W|}{n(n-1)} \sum_{i=1}^n \sum_{\substack{j=1 \\ j \neq i}}^n \mathbf{1}\{ d_{ij} \leq r \} \, e_{ij}(r) \]
dengan keterangan: - \(d_{ij}\) :
jarak teramati antara titik \(i\) dan
\(j\)
- \(r \geq 0\) : radius/jarak
tertentu
- \(e_{ij}(r)\) : faktor koreksi tepi
(edge correction weight)
- \(n\) : jumlah total titik
- \(|W|\) : luas jendela observasi
(area studi)
- \(\mathbf{1}\{d_{ij} \leq r\}\) :
fungsi indikator, bernilai 1 jika \(d_{ij}
\leq r\), dan 0 jika tidak
Intuisi
Dengan kata lain, \(\hat{K}(r)\) adalah ukuran jumlah titik tetangga rata-rata dalam jarak \(r\), setelah dikoreksi tepi dan distandardisasi.
# Hitung jarak tetangga terdekat
kplot <- nndist(nztrees)
# Histogram lebih menarik
hist(kplot,
breaks = 15, # jumlah kelas histogram
col = "skyblue", # warna batang
border = "white", # warna tepi batang
main = "Histogram Jarak Tetangga Terdekat (nztrees)",
xlab = "Jarak (satuan koordinat)",
ylab = "Frekuensi",
cex.main = 1, cex.lab = 0.9, cex.axis = 0.8)
# Tambahkan garis rata-rata
abline(v = mean(kplot), col = "red", lwd = 2, lty = 2)
text(mean(kplot), max(hist(kplot, plot = FALSE)$counts),
labels = paste0("Mean = ", round(mean(kplot), 2)),
pos = 4, col = "red", cex = 0.8)
# Tambahkan garis median
abline(v = median(kplot), col = "darkgreen", lwd = 2, lty = 3)
text(median(kplot), max(hist(kplot, plot = FALSE)$counts) * 0.9,
labels = paste0("Median = ", round(median(kplot), 2)),
pos = 4, col = "darkgreen", cex = 0.8)# Hitung K-function dengan koreksi Ripley
K <- Kest(nztrees, correction = "Ripley")
# Plot dengan judul
plot(K,
main = "Empirical K-Function (Ripley's Correction)",
xlab = "Jarak (r)",
ylab = "K(r)")## Generating 999 simulations of CSR ...
## 1, 2, 3, ......10.........20.........30.........40.........50.........60..
## .......70.........80.........90.........100.........110.........120.........130
## .........140.........150.........160.........170.........180.........190........
## .200.........210.........220.........230.........240.........250.........260......
## ...270.........280.........290.........300.........310.........320.........330....
## .....340.........350.........360.........370.........380.........390.........400..
## .......410.........420.........430.........440.........450.........460.........470
## .........480.........490.........500.........510.........520.........530........
## .540.........550.........560.........570.........580.........590.........600......
## ...610.........620.........630.........640.........650.........660.........670....
## .....680.........690.........700.........710.........720.........730.........740..
## .......750.........760.........770.........780.........790.........800.........810
## .........820.........830.........840.........850.........860.........870........
## .880.........890.........900.........910.........920.........930.........940......
## ...950.........960.........970.........980.........990........
## 999.
##
## Done.
# Plot dengan judul
plot(E,
main = "Envelope Empirical K-Function (nztrees)",
xlab = "Jarak (r)",
ylab = "K(r)")## Generating 999 simulations of CSR ...
## 1, 2, 3, ......10.........20.........30.........40.........50.........60..
## .......70.........80.........90.........100.........110.........120.........130
## .........140.........150.........160.........170.........180.........190........
## .200.........210.........220.........230.........240.........250.........260......
## ...270.........280.........290.........300.........310.........320.........330....
## .....340.........350.........360.........370.........380.........390.........400..
## .......410.........420.........430.........440.........450.........460.........470
## .........480.........490.........500.........510.........520.........530........
## .540.........550.........560.........570.........580.........590.........600......
## ...610.........620.........630.........640.........650.........660.........670....
## .....680.........690.........700.........710.........720.........730.........740..
## .......750.........760.........770.........780.........790.........800.........810
## .........820.........830.........840.........850.........860.........870........
## .880.........890.........900.........910.........920.........930.........940......
## ...950.........960.........970.........980.........990........
## 999.
##
## Done.
##
## Maximum absolute deviation test of CSR
## Monte Carlo test based on 999 simulations
## Summary function: K(r)
## Reference function: theoretical
## Alternative: two.sided
## Interval of distance values: [0, 23.75] feet
## Test statistic: Maximum absolute deviation
## Deviation = observed minus theoretical
##
## data: nztrees
## mad = 119.43, rank = 383, p-value = 0.383
Interpretasi
Garis hitam = \(\hat{K}_{obs}(r)\), fungsi K empiris dari data.
Garis merah = \(K_{theo}(r)\), fungsi K teoretis di bawah asumsi CSR (Complete Spatial Randomness).
Daerah abu-abu = envelope dari 999 simulasi CSR.
Statistik uji: MAD (Maximum Absolute Deviation) = 119.43
p-value = 0.391
Karena p-value = 0.391 > 0.05, maka kita gagal menolak hipotesis nol (H0). → Artinya, tidak ada bukti signifikan bahwa pola titik menyimpang dari CSR. Secara visual, garis \(\hat{K}_{obs}(r)\) (hitam) hampir selalu berada di dalam envelope abu-abu, dan dekat dengan garis teoretis merah → Ini menunjukkan pola titik konsisten dengan sebaran acak. Jika garis hitam berada di atas envelope → indikasi cluster. Jika garis hitam berada di bawah envelope → indikasi uniform/regular. Pada grafik ini, tidak ada deviasi signifikan.
Berdasarkan analisis fungsi K dan uji MAD dengan 999 simulasi,
sebaran pohon dalam dataset nztrees tidak berbeda
signifikan dari pola acak (CSR). Dengan kata lain, distribusi
titik dapat dianggap acak.
Nearest Neighbor Analysis menganalisis jarak tetangga terdekat antar titik. Ia mendeteksi apakah pola: - Cluster (jarak aktual cenderung lebih kecil dari acak), - Uniform/regular (jarak aktual cenderung lebih besar dari acak), - Acak/CSR (Complete Spatial Randomness).
Dua keluaran utama: 1. Distribusi jarak tetangga
terdekat: vektor jarak untuk setiap titik. 2. Indeks
Clark–Evans \(R\): rasio
antara rataan jarak aktual vs rataan teoretis di bawah CSR.
\[
R = \frac{\bar{r}_{\text{obs}}}{\bar{r}_{\text{CSR}}}, \quad
\bar{r}_{\text{CSR}} = \frac{1}{2\sqrt{\lambda}}, \
\lambda=\frac{n}{|W|}
\] Interpretasi: \(R<1\)
(cluster), \(R\approx1\) (acak), \(R>1\) (regular).
## Planar point pattern: 86 points
## window: rectangle = [0, 153] x [0, 95] feet
Hitung Jarak Tetangga Terdekat & Ringkasan
# vektor jarak tetangga terdekat (satuan sesuai koordinat data)
r_nn <- nndist(pp)
# ringkasan deskriptif
summary(r_nn)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 2.000 3.162 6.243 6.856 10.198 17.029
# intensitas (lambda) dan luas jendela |W|
n <- pp$n
A <- area.owin(Window(pp))
lam <- n / A
rbar_obs <- mean(r_nn)
rbar_csr <- 1 / (2*sqrt(lam)) # ekspektasi mean NN di bawah CSR
c(n = n, area = A, lambda = lam,
mean_NN_obs = rbar_obs, mean_NN_CSR = rbar_csr)## n area lambda mean_NN_obs mean_NN_CSR
## 8.600000e+01 1.453500e+04 5.916753e-03 6.855581e+00 6.500224e+00
Clark–Evans 𝑅 R & Uji CSR
# Clark-Evans test (tersedia di spatstat.explore/spatstat lama)
# method default "D" (Donelly) dengan edge correction
res_ce <- clarkevans(pp, correction = "D") # menghasilkan R dan p-value
res_ce## [1] 1.004879
R<1 dan p-value < 0.05 → cluster signifikan
\(R\approx1\) dan p-value ≥ 0.05 → konsisten acak
R>1 dan p-value < 0.05 → regular/uniform signifikan
Fungsi G(r) (Nearest Neighbor Distribution Function)
Fungsi G(r)= peluang jarak tetangga terdekat ≤r. Dibandingkan dengan \(G_{\text{theo}}(r)\) di bawah CSR
G <- Gest(pp, correction = "rs") # "rs" = Kaplan–Meier style edge correction
plot(G, main = "Fungsi G(r) — Nearest Neighbor Distribution (nztrees)",
xlab = "Jarak (r)", ylab = "G(r)")Kurva empiris naik lebih cepat dari teoretis → cluster.
Kurva empiris lebih lambat dari teoretis → regular
Dekat teoretis → acak.
Envelope Monte Carlo untuk G(r)
## Generating 199 simulations of CSR ...
## 1, 2, 3, 4.6.8.10.12.14.16.18.20.22.24.26.28.30.32.34
## .36.38.40.42.44.46.48.50.52.54.56.58.60.62.64.66.68.70.72.74
## .76.78.80.82.84.86.88.90.92.94.96.98.100.102.104.106.108.110.112.114
## .116.118.120.122.124.126.128.130.132.134.136.138.140.142.144.146.148.150.152.154
## .156.158.160.162.164.166.168.170.172.174.176.178.180.182.184.186.188.190.192.194
## .196.198
## 199.
##
## Done.
Berdasarkan analisis jarak tetangga terdekat (fungsi G, envelope
Monte Carlo, dan indeks Clark–Evans), pola sebaran pohon pada dataset
nztrees tidak berbeda signifikan dari pola acak
(Complete Spatial Randomness). Artinya, distribusi titik
cenderung acak dan tidak menunjukkan pengelompokan atau
keteraturan yang jelas.
Definisi: \[ J(r) \;=\; \frac{1 - F(r)}{\,1 - G(r)\,} \] dengan:
- \(G(r)\): Nearest-Neighbor Distribution — peluang jarak tetangga terdekat \(\le r\).
- \(F(r)\): Empty-Space Distribution — peluang jarak dari titik acak di area ke titik terdekat \(\le r\).
Di bawah CSR (Complete Spatial Randomness): \(G(r)=F(r)\) sehingga \(J(r)=1\) untuk semua \(r\).
## Planar point pattern: 86 points
## window: rectangle = [0, 153] x [0, 95] feet
J <- Jest(pp, correction = "rs")
# Plot dengan garis teori (CSR = 1)
plot(J,
main = "J-Function (nztrees)",
xlab = "Jarak (r)", ylab = "J(r)")
abline(h = 1, lty = 2) # garis teoretis CSRset.seed(123)
EJ <- envelope(pp, Jest, nsim = 199, correction = "rs") # bisa naikkan nsim (mis. 999)## Generating 199 simulations of CSR ...
## 1, 2, 3, 4.6.8.10.12.14.16.18.20.22.24.26.28.30.32.34
## .36.38.40.42.44.46.48.50.52.54.56.58.60.62.64.66.68.70.72.74
## .76.78.80.82.84.86.88.90.92.94.96.98.100.102.104.106.108.110.112.114
## .116.118.120.122.124.126.128.130.132.134.136.138.140.142.144.146.148.150.152.154
## .156.158.160.162.164.166.168.170.172.174.176.178.180.182.184.186.188.190.192.194
## .196.198
## 199.
##
## Done.
plot(EJ,
main = "Envelope J-Function — Monte Carlo (nztrees)",
xlab = "Jarak (r)", ylab = "J(r)")
abline(h = 1, lty = 2)# Maximum Absolute Deviation test (two-sided) terhadap fungsi teoretis (CSR = 1)
set.seed(123)
madJ <- mad.test(pp, Jest, nsim = 499) # tambah nsim untuk ketelitian (mis. 999)## Generating 499 simulations of CSR ...
## 1, 2, 3, .5....10....15....20....25....30....35....40....45....50..
## ..55....60....65....70....75....80....85....90....95....100....105....110
## ....115....120....125....130....135....140....145....150....155....160....165...
## .170....175....180....185....190....195....200....205....210....215....220....225.
## ...230....235....240....245....250....255....260....265....270....275....280....
## 285....290....295....300....305....310....315....320....325....330....335....340..
## ..345....350....355....360....365....370....375....380....385....390....395....400
## ....405....410....415....420....425....430....435....440....445....450....455...
## .460....465....470....475....480....485....490....495...
## 499.
##
## Done.
## Warning in envelopeTest(X, ..., exponent = Inf, alternative = alternative, :
## Some function values were infinite, NA or NaN at distances r up to
## 24.8872094375391 feet. Consider specifying a shorter 'rinterval'
## Warning in envelopeTest(X, ..., exponent = Inf, alternative = alternative, :
## Some deviation values were Inf, NA or NaN
## Warning in envelopeTest(X, ..., exponent = Inf, alternative = alternative, :
## Some simulated deviations were Inf, NA or NaN
##
## Maximum absolute deviation test of CSR
## Monte Carlo test based on 499 simulations
## Summary function: J(r)
## Reference function: theoretical
## Alternative: two.sided
## Interval of distance values: [0, 24.8872094375391] feet
## Test statistic: Maximum absolute deviation
## Deviation = observed minus theoretical
##
## data: pp
## mad = 15.459, rank = 132, p-value = 0.264
## Generating 499 simulations of CSR ...
## 1, 2, 3, .5....10....15....20....25....30....35....40....45....50..
## ..55....60....65....70....75....80....85....90....95....100....105....110
## ....115....120....125....130....135....140....145....150....155....160....165...
## .170....175....180....185....190....195....200....205....210....215....220....225.
## ...230....235....240....245....250....255....260....265....270....275....280....
## 285....290....295....300....305....310....315....320....325....330....335....340..
## ..345....350....355....360....365....370....375....380....385....390....395....400
## ....405....410....415....420....425....430....435....440....445....450....455...
## .460....465....470....475....480....485....490....495...
## 499.
##
## Done.
## Warning in envelopeTest(X, ..., exponent = 2, alternative = alternative, : Some
## function values were infinite, NA or NaN at distances r up to 24.8872094375391
## feet. Consider specifying a shorter 'rinterval'
## Warning in envelopeTest(X, ..., exponent = 2, alternative = alternative, : Some
## deviation values were Inf, NA or NaN
## Warning in envelopeTest(X, ..., exponent = 2, alternative = alternative, : Some
## simulated deviations were Inf, NA or NaN
##
## Diggle-Cressie-Loosmore-Ford test of CSR
## Monte Carlo test based on 499 simulations
## Summary function: J(r)
## Reference function: theoretical
## Alternative: two.sided
## Interval of distance values: [0, 24.8872094375391] feet
## Test statistic: Integral of squared absolute deviation
## Deviation = observed minus theoretical
##
## data: pp
## u = 162.39, rank = 124, p-value = 0.248
Hasil Uji Statistik
Maximum Absolute Deviation (MAD) test
Statistik uji = 15.459
p-value = 0.264 (> 0.05) → Gagal menolak H₀
Diggle–Cressie–Loosmore–Ford (DCLF) test
Statistik uji = 162.39
p-value = 0.248 (> 0.05) → Gagal menolak H₀
Kesimpulan
H₀: Pola sebaran titik mengikuti CSR (Complete Spatial Randomness). Karena p-value keduanya > 0.05, tidak ada bukti signifikan untuk menolak H₀.
Dengan kata lain, fungsi J menunjukkan bahwa pola sebaran
pohon nztrees konsisten dengan pola acak. Tidak
ada indikasi kuat adanya pengelompokan (cluster) maupun
pola regular/uniform. Berdasarkan J-Function dan uji
Monte Carlo (MAD & DCLF), sebaran pohon pada dataset
nztrees tidak berbeda signifikan dari pola acak
(CSR).
| Function / Metode | Konsep Dasar | Input Data | Skala Analisis | Kelebihan | Keterbatasan | Interpretasi Pola |
|---|---|---|---|---|---|---|
| Nearest Neighbor Analysis (NNA) | Mengukur jarak tetangga terdekat tiap titik, dibandingkan dengan distribusi acak | Koordinat titik | Skala lokal (titik ke titik) | Mudah dihitung, cepat, intuitif (Clark–Evans R, G-function) | Hanya fokus jarak terdekat → bisa bias jika pola berbeda di jarak menengah/besar | R < 1 → cluster, R ≈ 1 → acak, R > 1 → uniform |
| Metode Kuadran (Quadrat Count / VMR) | Membagi area jadi sel sama besar, lalu menganalisis ragam vs rata-rata hitungan titik per sel | Jumlah titik per sel | Skala global kasar | Sederhana, cocok untuk data lapangan | Sensitif pada ukuran sel, tidak multi-skala | VMR < 1 → uniform, VMR ≈ 1 → acak, VMR > 1 → cluster |
| K-Function (Ripley’s K / Besag’s L) | Menghitung jumlah rata-rata tetangga dalam radius r, distandardisasi dengan intensitas, dikoreksi tepi | Koordinat titik + jendela observasi | Multi-skala (berbagai jarak r) | Menangkap pola pada semua skala jarak, formal, tersedia uji Monte Carlo | Lebih kompleks, butuh edge correction | L(r) ≈ 0 → acak, L(r) > 0 → cluster, L(r) < 0 → uniform |
| J-Function | Kombinasi F(r) (empty space) dan G(r) (nearest neighbor): \(J(r) = \frac{1-F(r)}{1-G(r)}\) | Koordinat titik + jendela | Multi-skala, stabil | Lebih robust daripada F atau G sendiri, membandingkan jarak acak vs tetangga | Butuh simulasi acak untuk envelope | J(r) ≈ 1 → acak, J(r) < 1 → cluster, J(r) > 1 → uniform |
Untuk menguji Kuadran (Quadrat/VMR), K–Function (Ripley/Besag’s L), dan J–Function dengan data bangkitan, pola terbaik adalah:
pola nol (null) = CSR/Poisson (Complete Spatial Randomness),
pola alternatif-1 (clustered),
pola alternatif-2 (regular/inhibitory).
Dengan begitu kamu bisa mengecek ukuran uji (type-I error) di bawah CSR dan power saat pola bukan CSR. Bangkitan memiliki karakteristik
Intensitas target: mis. λ ≈ 1000 / (56×38) ≈ 0.47 titik/m² (atau pilih jumlah titik tetap, mis. N≈1000).
Null – CSR: Poisson
(rpoispp).
rThomas) → tentukan skala cluster scale (≈1
m), rata² anak per induk mu (≈10), dan
kappa = λ/mu.rStrauss / rHardcore) → pilih radius inhibisi
R (≈0.6–1 m) dan parameter interaksi
(gamma<1 untuk Strauss).# ---- Paket ----
library(spatstat.geom)
library(spatstat.random)
set.seed(42)
# ---- Window & intensitas ----
win <- owin(xrange = c(0, 56), yrange = c(0, 38))
area <- area.owin(win)
N_target <- 1000
lambda <- N_target / area # ≈ 0.047 titik/m2 (diperbaiki)
# ==== NULL: CSR ====
pp_csr <- rpoispp(lambda = lambda, win = win)
# ==== CLUSTER: Thomas ====
mu <- 5 # Diperkecil untuk menghindari terlalu banyak titik
scale <- 1
kappa <- lambda / mu
pp_cluster <- rThomas(kappa = kappa, scale = scale, mu = mu, win = win)
# ==== REGULAR: Alternatif lebih cepat ====
# Opsi 1: rSSI (Simple Sequential Inhibition) - lebih cepat
R <- 1.5 # Jarak minimum antar titik
n_regular <- N_target # Target jumlah titik
pp_regular <- rSSI(r = R, n = n_regular, win = win)## Warning in rSSI(r = R, n = n_regular, win = win): Gave up after 1000 attempts
## with only 605 points placed out of 1000
adjust_pattern <- function(pp, target = N_target) {
n <- npoints(pp)
if (n > target * 1.5) {
# Terlalu banyak titik, ambil sampel acak
indices <- sample(1:n, target)
pp <- pp[indices]
} else if (n < target * 0.5) {
# Terlalu sedikit titik, buat ulang dengan parameter berbeda
pp <- rpoispp(lambda = lambda, win = win)
}
return(pp)
}
pp_csr <- adjust_pattern(pp_csr)
pp_cluster <- adjust_pattern(pp_cluster)
pp_regular <- adjust_pattern(pp_regular)
# ---- Plot tiga pola dengan style lebih cantik ----
par(mfrow=c(1,3), mar=c(2,2,3,1), bg="white")
cols <- c("darkblue", "darkgreen", "firebrick")
plot(pp_csr,
pch=19, cols=cols[1], cex=0.6,
main=paste("CSR\n(n =", npoints(pp_csr), ")"),
cex.main=1, axes=FALSE, ann=FALSE)
box(col="gray40")
plot(pp_cluster,
pch=19, cols=cols[2], cex=0.6,
main=paste("Cluster\n(n =", npoints(pp_cluster), ")"),
cex.main=1, axes=FALSE, ann=FALSE)
box(col="gray40")
plot(pp_regular,
pch=19, cols=cols[3], cex=0.6,
main=paste("Regular\n(n =", npoints(pp_regular), ")"),
cex.main=1, axes=FALSE, ann=FALSE)
box(col="gray40")## Jumlah titik:
## CSR: 1043
## Cluster: 966
## Regular: 605
qx <- 10; qy <- 10
qt <- quadratcount(pp_csr, nx = qx, ny = qy)
plot(pp_csr, main = "CSR + Grid Kuadran"); plot(qt, add=TRUE, cex=0.8, col="red")# Uji chi-square terhadap CSR homogen (intensitas konstan)
qt_test_csr <- quadrat.test(pp_csr, nx = qx, ny = qy)
qt_test_cluster <- quadrat.test(pp_cluster, nx = qx, ny = qy)
qt_test_regular <- quadrat.test(pp_regular, nx = qx, ny = qy)
qt_test_csr##
## Chi-squared test of CSR using quadrat counts
##
## data: pp_csr
## X2 = 84.996, df = 99, p-value = 0.3181
## alternative hypothesis: two.sided
##
## Quadrats: 10 by 10 grid of tiles
##
## Chi-squared test of CSR using quadrat counts
##
## data: pp_cluster
## X2 = 434.62, df = 99, p-value < 2.2e-16
## alternative hypothesis: two.sided
##
## Quadrats: 10 by 10 grid of tiles
##
## Chi-squared test of CSR using quadrat counts
##
## data: pp_regular
## X2 = 20.62, df = 99, p-value < 2.2e-16
## alternative hypothesis: two.sided
##
## Quadrats: 10 by 10 grid of tiles
# Catatan: VMR = Var(Counts) / Mean(Counts)
vmr <- function(pp, nx=qx, ny=qy){
tab <- as.vector(quadratcount(pp, nx=nx, ny=ny))
var(tab)/mean(tab)
}
c(VMR_CSR = vmr(pp_csr), VMR_Cluster = vmr(pp_cluster), VMR_Regular = vmr(pp_regular))## VMR_CSR VMR_Cluster VMR_Regular
## 0.8585471 4.3901123 0.2082812
K_csr <- Kest(pp_csr, correction="border")
L_csr <- Lest(pp_csr, correction="border")
# Envelope terhadap CSR (NULL) — band pembanding
envK_csr <- envelope(pp_csr, fun=Kest, nsim=199, correction="border", global=TRUE)## Generating 199 simulations of CSR ...
## 1, 2, 3, 4.6.8.10.12.14.16.18.20.22.24.26.28.30.32.34
## .36.38.40.42.44.46.48.50.52.54.56.58.60.62.64.66.68.70.72.74
## .76.78.80.82.84.86.88.90.92.94.96.98.100.102.104.106.108.110.112.114
## .116.118.120.122.124.126.128.130.132.134.136.138.140.142.144.146.148.150.152.154
## .156.158.160.162.164.166.168.170.172.174.176.178.180.182.184.186.188.190.192.194
## .196.198
## 199.
##
## Done.
## Generating 199 simulations of CSR ...
## 1, 2, 3, 4.6.8.10.12.14.16.18.20.22.24.26.28.30.32.34
## .36.38.40.42.44.46.48.50.52.54.56.58.60.62.64.66.68.70.72.74
## .76.78.80.82.84.86.88.90.92.94.96.98.100.102.104.106.108.110.112.114
## .116.118.120.122.124.126.128.130.132.134.136.138.140.142.144.146.148.150.152.154
## .156.158.160.162.164.166.168.170.172.174.176.178.180.182.184.186.188.190.192.194
## .196.198
## 199.
##
## Done.
# Untuk data alternatif, bandingkan ke envelope CSR:
envK_cluster <- envelope(pp_cluster, fun=Kest, nsim=199, simulate=expression(rpoispp(lambda, win=win)),
correction="border", global=TRUE)## Generating 398 simulations by evaluating expression (199 to estimate the mean
## and 199 to calculate envelopes) ...
## 1, 2, 3, .5....10....15....20....25....30....35....40....45....50..
## ..55....60....65....70....75....80....85....90....95....100....105....110
## ....115....120....125....130....135....140....145....150....155....160....165...
## .170....175....180....185....190....195....200....205....210....215....220....225.
## ...230....235....240....245....250....255....260....265....270....275....280....
## 285....290....295....300....305....310....315....320....325....330....335....340..
## ..345....350....355....360....365....370....375....380....385....390....395..
## 398.
##
## Done.
envK_regular <- envelope(pp_regular, fun=Kest, nsim=199, simulate=expression(rpoispp(lambda, win=win)),
correction="border", global=TRUE)## Generating 398 simulations by evaluating expression (199 to estimate the mean
## and 199 to calculate envelopes) ...
## 1, 2, 3, .5....10....15....20....25....30....35....40....45....50..
## ..55....60....65....70....75....80....85....90....95....100....105....110
## ....115....120....125....130....135....140....145....150....155....160....165...
## .170....175....180....185....190....195....200....205....210....215....220....225.
## ...230....235....240....245....250....255....260....265....270....275....280....
## 285....290....295....300....305....310....315....320....325....330....335....340..
## ..345....350....355....360....365....370....375....380....385....390....395..
## 398.
##
## Done.
par(mfrow=c(1,3), mar=c(4,4,2,1))
plot(envK_csr, main="Envelope K – data = CSR")
plot(envK_cluster, main="Envelope K – data = Cluster")
plot(envK_regular, main="Envelope K – data = Regular")# Di bawah CSR, J(r) ≈ 1. Di atas 1 → regular; di bawah 1 → cluster.
J_csr <- Jest(pp_csr)
J_cluster <- Jest(pp_cluster)
J_regular <- Jest(pp_regular)
par(mfrow=c(1,3), mar=c(4,4,2,1))
plot(J_csr, main="J – CSR"); abline(h=1, lty=2)
plot(J_cluster, main="J – Cluster"); abline(h=1, lty=2)
plot(J_regular, main="J – Regular"); abline(h=1, lty=2)nn_summary <- function(pp, name){
lam_hat <- intensity(pp)
rbar_obs <- mean(nndist(pp))
rbar_csr <- 1/(2*sqrt(lam_hat))
R_ce <- rbar_obs / rbar_csr
ce <- clarkevans.test(pp, correction = "Donnelly") # uji dengan koreksi tepi
data.frame(
Pattern = name,
n = npoints(pp),
lambda_hat = lam_hat,
mean_NND = rbar_obs,
expected_NND_CSR = rbar_csr,
R_ClarkEvans = R_ce,
CE_p_value = ce$p.value,
CE_alternative = ce$alternative
)
}
tab_nna <- rbind(
nn_summary(pp_csr, "CSR (Poisson)"),
nn_summary(pp_cluster, "Cluster (Thomas)"),
nn_summary(pp_regular, "Regular (Strauss)")
)
print(tab_nna, row.names = FALSE)## Pattern n lambda_hat mean_NND expected_NND_CSR R_ClarkEvans
## CSR (Poisson) 1043 0.4901316 0.7276131 0.7141898 1.0187951
## Cluster (Thomas) 966 0.4539474 0.6349312 0.7421082 0.8555776
## Regular (Strauss) 605 0.2843045 1.6302525 0.9377307 1.7385082
## CE_p_value CE_alternative
## 0.7339719 two-sided
## 0.0000000 two-sided
## 0.0000000 two-sided
# ===== Visualisasi distribusi jarak NN (opsional, biar intuitif) =====
par(mfrow=c(1,3), mar=c(4,4,2,1))
hist(nndist(pp_csr), breaks=30, main="Nearest-Neighbor Distance\nCSR", xlab="distance", col="gray90", border="white")
abline(v = 1/(2*sqrt(intensity(pp_csr))), lty=2)
hist(nndist(pp_cluster), breaks=30, main="Nearest-Neighbor Distance\nCluster", xlab="distance", col="gray90", border="white")
abline(v = 1/(2*sqrt(intensity(pp_cluster))), lty=2)
hist(nndist(pp_regular), breaks=30, main="Nearest-Neighbor Distance\nRegular", xlab="distance", col="gray90", border="white")
abline(v = 1/(2*sqrt(intensity(pp_regular))), lty=2)CSR (Poisson): p-value = 0.318 → tidak
menolak CSR.
VMR ≈ 0.86 mendekati 1 → sebaran acak.
Cluster (Thomas): p-value < 2.2e-16 →
tolak CSR.
VMR ≈ 4.39 >> 1 → indikasi kuat adanya pengelompokan.
Regular (Strauss): p-value < 2.2e-16 →
tolak CSR.
VMR ≈ 0.21 << 1 → titik lebih merata dari CSR, sesuai pola
inhibisi.
Hasil dari keempat metode konsisten dan saling melengkapi. Dataset CSR (Poisson) menunjukkan distribusi acak baik secara global (VMR ~1, K/L dalam envelope) maupun lokal (J ~1, R ~1). Dataset Thomas process (Cluster) jelas memperlihatkan pola pengelompokan, ditunjukkan oleh VMR tinggi, K-Function di atas envelope, J < 1, dan nilai R < 1 yang signifikan. Sebaliknya, dataset Strauss process (Regular) menegaskan adanya inhibisi antar titik, tercermin dari VMR < 1, K-Function di bawah envelope, J > 1, dan R jauh lebih besar dari 1. Dengan demikian, kombinasi keempat pendekatan ini efektif untuk membedakan pola spasial acak, mengelompok, dan teratur, sekaligus menunjukkan bahwa uji-uji ini saling melengkapi: Quadrat lebih global, K/L multiskala, J berfokus pada tetangga terdekat, dan NNA memberi ukuran ringkas dengan uji signifikansi yang jelas.
Sumber
Diggle, P. J. (2013). Statistical Analysis of Spatial and Spatio-Temporal Point Patterns (3rd ed.). CRC Press.
Ripley, B. D. (1976). The second-order analysis of stationary point processes. Journal of Applied Probability, 13(2), 255–266.
Besag, J. (1977). Contribution to the discussion of Dr. Ripley’s paper. Journal of the Royal Statistical Society. Series B (Methodological), 39(2), 193–195.
van Lieshout, M. N. M., & Baddeley, A. J. (1996). A nonparametric measure of spatial interaction in point patterns. Statistica Neerlandica, 50(3), 344–361.
Clark, P. J., & Evans, F. C. (1954). Distance to nearest neighbor as a measure of spatial relationships in populations. Ecology, 35(4), 445–453.
Cressie, N. A. C. (1993). Statistics for Spatial Data. Wiley.