Analisis Perbandingan Kuadran, K-Function, J-Function, dan NNA dalam Identifikasi Pola Spasial

I Gusti Ngurah Sentana Putra (M0501241019)

Tugas

Pada tugas ini, diminta untuk:

  1. Melakukan analisis terhadap data titik (boleh dari R package atau sumber lain)

  2. Menjelaskan hasil analisis (minimal mencakup metode kuadran dan fungsi-K)

  3. Menuliskannya dan mengumpulkan dalam bentuk html, lengkap dengan keterangannya.

Data Package R

Data

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.

library(spatstat)
## 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'
data(nztrees)
nztrees
## Planar point pattern: 86 points
## window: rectangle = [0, 153] x [0, 95] feet
summary(nztrees)
## 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
plot(nztrees, main = "Lokasi Pohon Beech di New Zealand")

nztrees <- nztrees
contour(density(nztrees,5),axes = F)

Metode Kuadran

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:

  1. Bagi area studi menjadi \(m\) sel dengan ukuran sama (ditentukan sesuai skala).
  2. Hitung rata-rata jumlah titik per sel: \[ \bar{x} = \frac{N}{m} \] dengan \(N =\) total banyak titik, \(m =\) jumlah sel.
  3. Hitung ragam jumlah titik per sel: \[ S^2 = \frac{\sum_{i=1}^m (x_i - \bar{x})^2}{m-1} \] dengan \(x_i =\) banyak titik pada sel ke-\(i\).
  4. Hitung rasio ragam terhadap rata-rata (Variance to Mean Ratio / VMR): \[ VMR = \frac{S^2}{\bar{x}} \]

Interpretasi VMR

  • \(VMR \approx 0\) : titik tersebar uniform (sistematis).
  • \(VMR \approx 1\) : titik tersebar acak.
  • \(VMR > 1\) : titik tersebar mengelompok (cluster).

Hipotesis

  • \(H_0\): Pola sebaran titik acak.
  • \(H_1\): Pola sebaran titik tidak acak.

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\%\))

  • Jika \(Z > 1.96\) : Tolak \(H_0\), kesimpulan cluster.
  • Jika \(Z < -1.96\) : Tolak \(H_0\), kesimpulan uniform.
  • Jika \(-1.96 \leq Z \leq 1.96\) : Terima \(H_0\), kesimpulan acak.

Aplikasi dalam R

q <- quadratcount(nztrees, nx = 4, ny = 3)
q
##              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
plot(q)

mu <- mean(q)
sigma <- sd(q)^2
VMR <- sigma/mu
VMR
## [1] 0.807611
quadrat.test(q)
## 
##  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
quadrat.test(q, alt = "regular")
## 
##  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.

Empirical K-Function

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

  • \(\hat{K}(r)\) menyatakan rata-rata kumulatif banyaknya titik yang berada dalam jarak \(r\) dari sebuah titik acak.
  • Fungsi ini dikoreksi terhadap efek tepi (agar titik dekat batas area tidak bias).
  • Nilai juga distandardisasi dengan intensitas (jumlah titik per unit area).

Dengan kata lain, \(\hat{K}(r)\) adalah ukuran jumlah titik tetangga rata-rata dalam jarak \(r\), setelah dikoreksi tepi dan distandardisasi.

Aplikasi Dalam R

kplot <- nndist(nztrees)
hist(kplot)

# 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)")

E<-envelope(nztrees,Kest, nsim=999)
## 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)")

mad.test(nztrees, Kest, nsim=999, alternative="two.sided")
## 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 (NNA)

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).

Aplikasi dalam R

data(nztrees)   # objek 'ppp'
pp <- nztrees
pp
## 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)

set.seed(123)
envG <- envelope(pp, Gest, nsim = 199, correction = "rs")  # 199 simulasi CSR
## 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(envG, main = "Envelope G(r) — Monte Carlo (nztrees)",
     xlab = "Jarak (r)", ylab = "G(r)")

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.

J-Function

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\).

Aplikasi Dalam R

data(nztrees)
pp <- nztrees
pp
## 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 CSR

set.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)

  • Jika kurva \(\hat J(r)\) < 1 di luar envelope → indikasi cluster.
  • Jika kurva \(\hat J(r)\) > 1 di luar envelope → indikasi regular/uniform.
  • Jika tetap di dalam envelope → konsisten acak (CSR).
# 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
madJ
## 
##  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
set.seed(123)
dclfJ <- dclf.test(pp, Jest, nsim = 499)
## 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
dclfJ
## 
##  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

  1. Maximum Absolute Deviation (MAD) test

    • Statistik uji = 15.459

    • p-value = 0.264 (> 0.05) → Gagal menolak H₀

  2. 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).

Perbandingan Metode Analisis Pola Titik

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

Kajian Simulasi Perbandingan Metode

Rancangan Simulasi

Untuk menguji Kuadran (Quadrat/VMR), K–Function (Ripley/Besag’s L), dan J–Function dengan data bangkitan, pola terbaik adalah:

  1. pola nol (null) = CSR/Poisson (Complete Spatial Randomness),

  2. pola alternatif-1 (clustered),

  3. 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).

  • Alternatif (Cluster): Thomas process (rThomas) → tentukan skala cluster scale (≈1 m), rata² anak per induk mu (≈10), dan kappa = λ/mu.
  • Alternatif (Regular): Strauss / Hardcore (rStrauss / rHardcore) → pilih radius inhibisi R (≈0.6–1 m) dan parameter interaksi (gamma<1 untuk Strauss).

Data

# ---- 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")

par(mfrow=c(1,1))


# ---- Ringkasan ----
cat("Jumlah titik:\n")
## Jumlah titik:
cat("CSR:", npoints(pp_csr), "\n")
## CSR: 1043
cat("Cluster:", npoints(pp_cluster), "\n")
## Cluster: 966
cat("Regular:", npoints(pp_regular), "\n")
## Regular: 605

Metode QUADRAT

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
qt_test_cluster
## 
##  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
qt_test_regular
## 
##  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-Function & Besag’s L

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.
envL_csr <- envelope(pp_csr, fun=Lest, 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.
# 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")

par(mfrow=c(1,1))

J Function

# 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)

par(mfrow=c(1,1))

NNA

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)

par(mfrow=c(1,1))

Interpretasi

1. Metode Kuadran (Chi-Square & VMR)

  • 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.

2. K-Function (Ripley’s K / Besag’s L)

  • CSR: kurva observasi berada dalam envelope → konsisten dengan acak.
  • Cluster: K_obs di atas envelope → kelebihan pasangan titik pada jarak tertentu → indikasi cluster.
  • Regular: K_obs di bawah envelope → kekurangan pasangan titik → pola inhibisi/regular.

3. J-Function

  • CSR: J ≈ 1 di hampir semua jarak → sesuai acak.
  • Cluster: J < 1 → peluang tetangga dekat lebih besar → mendukung pola cluster.
  • Regular: J > 1 → peluang tetangga dekat lebih kecil → mendukung pola regular/inhibisi.

4. Nearest Neighbor Analysis (Clark–Evans R)

  • CSR: R ≈ 1.02, p-value 0.73 → sesuai CSR
  • Cluster: R ≈ 0.86, p-value < 0.001 → signifikan cluster.
  • Regular: R ≈ 1.74, p-value < 0.001 → signifikan regular.

Kesimpulan Simulasi

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.