DAPH Guide DAPH Guide

1บทนำและความสำคัญทางสาธารณสุขIntroduction and Public-Health Relevance

Spatial regression ขยาย linear regression ปกติให้สามารถจัดการกับสมมุติฐานที่ถูกละเมิดบ่อยที่สุดในข้อมูลภูมิศาสตร์ — นั่นคือ spatial autocorrelation (พื้นที่ที่ใกล้กันมักมีค่าคล้ายกัน) ซึ่งทำให้ OLS มี standard errors ผิด, t-statistics ที่บอกความนัยสำคัญเกินจริง, และค่าสัมประสิทธิ์ที่อาจ bias. DAPH รองรับ 3 รูปแบบหลัก: OLS (baseline + Moran's I และ LM diagnostics), Spatial Lag Model (SLM) ที่เพิ่ม spatial lag ของ outcome (ρWy) และ Spatial Error Model (SEM) ที่ปล่อย error term ให้มี autocorrelation (u = λWu + ε).

Spatial regression extends ordinary linear regression to handle the most-frequently-violated assumption in geographic data — spatial autocorrelation (nearby places tend to have similar values). When ignored, OLS yields wrong standard errors, inflated t-statistics, and potentially biased coefficients. DAPH supports the three canonical specifications: OLS (baseline + Moran's I and LM diagnostics), Spatial Lag Model (SLM) which adds a spatial lag of the outcome (ρWy), and Spatial Error Model (SEM) which lets the error term carry the spatial structure (u = λWu + ε).

2หลักการและสูตรLogic and Formula

กำหนดให้ y = outcome vector (n×1), X = design matrix (n×k), และ W = spatial weights matrix ที่ row-standardised. ทั้ง 3 แบบจำลองมีรูปแบบดังนี้:

Let y = outcome vector (n×1), X = design matrix (n×k), and W = row-standardised spatial weights matrix. The three specifications are:

Spatial Pseudo R² ของ DAPH ใช้สูตรมาตรฐานสำหรับทั้ง SLM และ SEM ตามที่ Anselin (1988) นิยามไว้:

DAPH's Spatial Pseudo R² uses the standard formula for both SLM and SEM as defined by Anselin (1988):

3สคริปต์ R และการทวนสอบกับ DAPHR Script and Verification against DAPH

ใช้ Columbus dataset (ตัวอย่าง canonical จาก Anselin 1988, มากับ R package spdep) — outcome คืออัตราอาชญากรรม (CRIME), predictors คือ HOVAL (มูลค่าบ้าน) และ INC (รายได้):

Using the canonical Columbus dataset (Anselin 1988, bundled with R package spdep) — outcome is CRIME rate, predictors HOVAL (housing value) and INC (income):

library(spdep)
library(spatialreg)

# Built-in Columbus dataset (49 neighbourhoods, Anselin 1988)
data(columbus, package = "spData")
col.gal.nb <- read.gal(system.file("etc/weights/columbus.gal", package = "spData"))
listw <- nb2listw(col.gal.nb, style = "W")   # row-standardised

# 1. OLS baseline (for comparison)
ols <- lm(CRIME ~ HOVAL + INC, data = columbus)
summary(ols)
lm.morantest(ols, listw)            # Moran's I on residuals
lm.LMtests(ols, listw, test = "all") # LM diagnostics

# 2. Spatial Lag Model (SLM)
slm <- lagsarlm(CRIME ~ HOVAL + INC, data = columbus, listw = listw)
summary(slm)

# 3. Spatial Error Model (SEM)
sem <- errorsarlm(CRIME ~ HOVAL + INC, data = columbus, listw = listw)
summary(sem)

# OLS:
#   (Intercept)  68.6190   HOVAL  -0.2739*   INC  -1.5973***
#   R² = 0.5524    Moran's I on residuals = 0.2226, p = 0.0021

# SLM:
#   rho (ρ)      = 0.4233     LR test p = 0.0035
#   HOVAL        = -0.2656*   INC = -1.0316**
#   AIC          = 374.8

# SEM:
#   lambda (λ)   = 0.5604     LR test p = 0.0011
#   HOVAL        = -0.2929*   INC = -0.9405*
#   AIC          = 374.1

* p<.05, ** p<.01, *** p<.001 — ค่าข้างต้นเป็นค่ามาตรฐานที่อ้างใน Anselin (1988) บทที่ 12 และ Bivand & Wong (2018).* p<.05, ** p<.01, *** p<.001 — values are the canonical figures reported in Anselin (1988) Ch.12 and Bivand & Wong (2018).

ค่า Columbus มาจาก Anselin (1988) Ch.12 (ทวนสอบกับ R spatialreg 1.3-x ตามที่ Bivand & Wong 2018 รายงาน) · ค่า TB-2021 Pseudo R² ทวนสอบบน dataset จริง.Columbus values from Anselin (1988) Ch.12 (re-verified against R spatialreg 1.3-x as reported in Bivand & Wong 2018). TB-2021 Pseudo R² values verified on the real dataset.

4การแปลผลและตัวอย่างInterpretation and Example

ขั้นที่ 1 — ตรวจว่าควรใช้ spatial model หรือไม่: ดู Moran's I บน OLS residuals — ถ้า p < .05 แปลว่า OLS ไม่พอ. จากนั้นใช้ LM diagnostics: ถ้า LM-lag มีนัยสำคัญและ LM-error ไม่มี → เลือก SLM; ตรงข้าม → เลือก SEM; ถ้ามีนัยสำคัญทั้งคู่ → เปรียบเทียบ robust LM (LM-lag-r vs LM-error-r) แล้วเลือกตัวที่ยังมีนัยสำคัญ.

Step 1 — Decide whether spatial is needed: check Moran's I on OLS residuals — if p < .05, OLS is inadequate. Then use LM diagnostics: significant LM-lag and non-significant LM-error → choose SLM; the opposite → SEM; both significant → compare robust LMs (LM-lag-r vs LM-error-r), pick whichever stays significant.

ขั้นที่ 2 — แปลค่า ρ และ λ: ทั้งสองตัวอยู่ในช่วง −1 ถึง 1. ρ > 0 ใน SLM แปลว่ามี spillover — ค่าใน outcome กระจายไปยังพื้นที่เพื่อนบ้าน (เช่น TB ในอำเภอเพื่อนบ้านสูง → TB ในอำเภอเราก็มีแนวโน้มสูง). λ > 0 ใน SEM แปลว่าปัจจัยที่ไม่ได้วัดมีโครงสร้างเชิงพื้นที่ — ทำให้ error correlate กันในพื้นที่ใกล้ๆ.

Step 2 — Interpret ρ and λ: both lie in (−1, 1). ρ > 0 in SLM means spillover — outcome values diffuse to neighbours (high TB in neighbouring districts → higher TB locally). λ > 0 in SEM means unobserved factors carry a spatial structure that correlates errors across nearby units.

ขั้นที่ 3 — เปรียบเทียบ β: ใน SLM/SEM ค่า β มักเล็กกว่า OLS เพราะ "การกระจายเชิงพื้นที่" ถูกแยกออกไปแล้ว · ใน SLM ต้องระวัง: ผลของ X ต่อ y ไม่ใช่ β เท่านั้น แต่เป็น (I − ρW)⁻¹ β ซึ่งรวม direct + indirect effect ผ่าน feedback loop.

Step 3 — Compare βs: SLM/SEM betas are usually smaller than OLS because spatial variance is partitioned out. In SLM, the impact of X on y is NOT just β — it is (I − ρW)⁻¹ β, the sum of direct + indirect effects through the spatial feedback.

1บทนำและความสำคัญทางสาธารณสุขIntroduction and Public-Health Relevance

Spatial Durbin Model (SDM) ขยาย Spatial Lag Model โดยเพิ่ม spatial lag ของตัวแปรอธิบาย (W·X) เข้าไปในสมการ — outcome ของพื้นที่หนึ่งจึงขึ้นกับทั้ง outcome ของเพื่อนบ้าน (ρWy) และ "ลักษณะของเพื่อนบ้าน" (WXθ) พร้อมกัน. รูปแบบนี้เป็นแบบจำลองที่ครอบคลุมที่สุดในตระกูล lag/error เพราะทั้ง SLM และ SEM เป็นกรณีพิเศษของ SDM. ใน DAPH, SDM เป็น tickbox ในหน้าเลือกแบบจำลองของ Spatial Regression (OLS / SLM / SEM / SDM) และรายงานจะแสดงบล็อกสัมประสิทธิ์ W×X พร้อมตาราง effects (direct / spillover / total).

The Spatial Durbin Model (SDM) extends the Spatial Lag Model by adding spatially lagged predictors (W·X) — an area's outcome depends simultaneously on its neighbours' outcome (ρWy) and on its neighbours' characteristics (WXθ). It is the most general member of the lag/error family: both SLM and SEM are special cases of SDM. In DAPH, SDM is a tickbox in the Spatial Regression model picker (OLS / SLM / SEM / SDM), and the report shows the W×X coefficient block plus the effects table (direct / spillover / total).

2หลักการและสูตรLogic and Formula

ช่วงความเชื่อมั่น 95% ของ direct / spillover / total ได้จากการจำลอง Monte-Carlo 1,000 รอบ ของพารามิเตอร์ที่ประมาณได้ (ρ, β, θ) ตามวิธีของ LeSage & Pace (2009).

95% CIs for direct / spillover / total come from 1,000 Monte-Carlo simulations of the fitted parameters (ρ, β, θ), following LeSage & Pace (2009).

3สคริปต์ R และการทวนสอบกับ DAPHR Script and Verification against DAPH

ทวนสอบบน TB dataset ราย 77 จังหวัด (queen contiguity, row-standardised W) — outcome คืออัตราป่วย TB, predictors คือ Workload และ Pop:

Verified on the 77-district TB dataset (queen contiguity, row-standardised W) — outcome TB rate, predictors Workload and Pop:

library(sf); library(spdep); library(spatialreg)
gdf <- st_read("districts.shp")
nb  <- poly2nb(gdf, queen = TRUE)
lw  <- nb2listw(nb, style = "W", zero.policy = TRUE)
m   <- lagsarlm(TB ~ Workload + Pop, data = gdf, listw = lw,
                Durbin = TRUE, zero.policy = TRUE)
summary(m)
impacts(m, listw = lw)

LeSage–Pace effects ทั้งสองฝั่งคำนวณจาก (I − ρW)⁻¹ เต็มเมทริกซ์ — ตรงกันถึงทศนิยม 8 ตำแหน่ง.LeSage–Pace effects on both sides computed from the full (I − ρW)⁻¹ matrix — identical to 8 decimal places.

4การแปลผลและตัวอย่างInterpretation and Example

ขั้นที่ 1 — ทดสอบว่า SDM จำเป็นหรือไม่: ใช้ Wald test (θ = 0) — ถ้าไม่ปฏิเสธ ให้ลดรูปเป็น SLM; ใช้ common-factor test (θ + ρβ = 0) — ถ้าไม่ปฏิเสธ ให้ลดรูปเป็น SEM. ถ้าปฏิเสธทั้งสอง ⇒ คง SDM แล้วรายงานผ่านตาราง effects เท่านั้น.

Step 1 — Test whether SDM is needed: the Wald test (θ = 0) — if not rejected, reduce to SLM; the common-factor test (θ + ρβ = 0) — if not rejected, reduce to SEM. If both are rejected ⇒ keep SDM and report through the effects table only.

ขั้นที่ 2 — อ่านตาราง effects ไม่ใช่ β ดิบ: ใน SDM ห้ามแปล β หรือ θ ตรง ๆ เพราะผลของตัวแปรหนึ่งไหลผ่าน feedback (I − ρW)⁻¹ ทั้งระบบ. Direct = ผลภายในพื้นที่ (รวม feedback), Spillover = ผลที่ส่งผ่านเพื่อนบ้าน, Total = ผลรวม — พร้อม 95% CI จาก Monte-Carlo 1,000 รอบ.

Step 2 — Read the effects table, not the raw βs: in SDM, never interpret β or θ directly — a predictor's effect propagates through the (I − ρW)⁻¹ feedback of the whole system. Direct = within-area effect (with feedback), Spillover = effect transmitted through neighbours, Total = their sum — each with a 95% CI from 1,000 Monte-Carlo simulations.

1บทนำและความสำคัญทางสาธารณสุขIntroduction and Public-Health Relevance

Spatial Panel Regression วิเคราะห์ข้อมูล พื้นที่ × เวลา ในรูปแบบ long format (1 แถวต่อพื้นที่ต่อช่วงเวลา, ต้องเป็น balanced panel — ทุกพื้นที่มีครบทุกช่วงเวลา). DAPH รองรับ 5 แบบจำลอง: Pooled OLS (baseline), FE-Lag, RE-Lag, FE-Error, RE-Error. Fixed effects (FE) ดูดซับคุณลักษณะประจำพื้นที่ที่ไม่เปลี่ยนตามเวลา (ใช้ within-area variation — จึงประมาณค่า predictor ที่คงที่ตามเวลาไม่ได้), ส่วน Random effects (RE) ใช้ทั้ง within + between variation และรายงาน φ (สัดส่วน variance ที่อยู่ระดับพื้นที่).

Spatial Panel Regression analyses areas × time data in long format (one row per area per period; a balanced panel is required — every area observed in every period). DAPH fits five models: Pooled OLS (baseline), FE-Lag, RE-Lag, FE-Error, RE-Error. Fixed effects (FE) absorb time-invariant area traits (within-area variation only — time-invariant predictors cannot be estimated), while Random effects (RE) use within + between variation and report φ (the share of variance at the area level).

2หลักการและสูตรLogic and Formula

การเลือกแบบจำลองใน DAPH: Hausman-type test ชี้ FE หรือ RE ก่อน จากนั้น DAPH เลือกแบบจำลองที่ใช้แปลผล (★) เป็นชนิดที่ Hausman แนะนำและมี AIC ต่ำสุดในกลุ่มนั้น.

Model selection in DAPH: the Hausman-type test first points to FE or RE, then DAPH marks as the interpreted model (★) the Hausman-preferred type with the lowest AIC.

3สคริปต์ R และการทวนสอบกับ DAPHR Script and Verification against DAPH

ทวนสอบบน balanced panel จำลอง 25 พื้นที่ × 4 ช่วงเวลา (queen weights บนพื้นที่ unique):

Verified on a synthetic balanced panel — 25 areas × 4 periods (queen weights over the unique areas):

library(splm); library(spdep)
# lw = queen listw over the UNIQUE areas
fe_lag <- spml(y ~ x1 + x2, data = pdata, index = c("area", "year"),
               listw = lw, model = "within", lag = TRUE,
               spatial.error = "none")
re_lag <- spml(y ~ x1 + x2, data = pdata, index = c("area", "year"),
               listw = lw, model = "random", lag = TRUE,
               spatial.error = "none")
sphtest(y ~ x1 + x2, data = pdata, index = c("area", "year"),
        listw = lw, spatial.model = "lag", method = "ML")

SE ของ RE-Lag ต่างกันเฉพาะธรรมเนียมการประมาณ variance components (≤ 0.0004 บน slopes) — point estimates ตรงกันทุกหลัก · Hausman ของ DAPH ใช้ quadratic form FE-vs-RE ซึ่งได้ค่าติดลบบนชุดข้อมูลนี้ = "ไม่มีหลักฐานคัดค้าน RE" (ธรรมเนียมมาตรฐาน) สอดคล้องกับ sphtest ของ R ที่สรุป RE เช่นกัน.RE-Lag SEs differ only in the variance-components convention (≤ 0.0004 on slopes) — point estimates identical to every digit. DAPH's FE-vs-RE quadratic form yields a negative statistic on this dataset = "no evidence against RE" (standard convention), agreeing with R's sphtest, which also concludes RE.

รอบแรกใช้ข้อมูลจำลองที่ไม่มี spatial lag จริง จึงยังไม่ได้ทดสอบว่า ρ ที่ประมาณได้ตรงกับ R หรือไม่ และยังไม่มีตัวแปรที่คงที่ตามเวลา. รอบนี้จึงจำลองใหม่ให้ครบเงื่อนไข: 100 พื้นที่ (ตาราง 10×10, Queen contiguity) × 8 ช่วงเวลา = 800 แถว จากสมการจริง y_it = a_i + 0.60·x1 − 0.35·x2 + 0.40·z_i + 0.35·W·y + e โดย z คงที่ภายในพื้นที่. ทั้งสองโปรแกรมได้รับ เมทริกซ์ W ชุดเดียวกัน และตารางข้อมูลชุดเดียวกัน.

The first round used data simulated without a real spatial lag, so it never tested whether the estimated ρ matches R, and it contained no time-invariant predictor. This round covers both: 100 areas (10×10 grid, Queen contiguity) × 8 periods = 800 rows generated from y_it = a_i + 0.60·x1 − 0.35·x2 + 0.40·z_i + 0.35·W·y + e, where z is constant within each area. Both programs received the same W matrix and the same data table.

ผลต่างสูงสุดที่วัดได้ทั้งตาราง = 1.3 × 10⁻⁴ (λ ของ FE-Error) ส่วนที่เหลืออยู่ระดับ 10⁻⁶–10⁻⁸ ซึ่งเป็นความต่างของอัลกอริทึมหาค่าเหมาะที่สุด ไม่ใช่ความต่างของสูตร.The largest measured difference across the whole table is 1.3 × 10⁻⁴ (FE-Error λ); everything else sits at 10⁻⁶–10⁻⁸, which is optimiser tolerance rather than a formula difference.

4การแปลผลและตัวอย่างInterpretation and Example

FE หรือ RE? ใช้ Hausman-type test บนบล็อก β ที่ใช้ร่วมกัน — p < .05 ⇒ RE ไม่ consistent ให้ใช้ FE; p ≥ .05 หรือสถิติติดลบ ⇒ ไม่มีหลักฐานคัดค้าน RE ซึ่งมี efficiency สูงกว่าและประมาณ predictor คงที่ตามเวลาได้. Lag หรือ Error? เทียบ AIC ภายในชนิดที่ Hausman แนะนำ — DAPH ทำให้อัตโนมัติและติดดาว ★ แบบจำลองที่ใช้แปลผล.

FE or RE? Use the Hausman-type test on the shared β block — p < .05 ⇒ RE is inconsistent, use FE; p ≥ .05 or a negative statistic ⇒ no evidence against RE, which is more efficient and can estimate time-invariant predictors. Lag or Error? Compare AIC within the Hausman-preferred type — DAPH automates this and stars ★ the interpreted model.

การอ่านค่า: ใน FE, β คือผลของการเปลี่ยนแปลง ภายในพื้นที่เดียวกัน ตามเวลา (ควบคุมคุณลักษณะถาวรของพื้นที่ทั้งหมดแล้ว). ใน RE, รายงาน φ ควบคู่ — φ สูงแปลว่า variance ส่วนใหญ่อยู่ระดับพื้นที่. ρ (Lag) แปลเป็น spillover ของ outcome ระหว่างพื้นที่ในช่วงเวลาเดียวกัน; λ (Error) คือโครงสร้างเชิงพื้นที่ของปัจจัยที่ไม่ได้วัด.

Reading the estimates: under FE, β is the effect of within-area change over time (all permanent area traits controlled). Under RE, report φ alongside — a high φ means most variance sits at the area level. ρ (Lag) is same-period outcome spillover between areas; λ (Error) is the spatial structure of unmeasured factors.

1บทนำและความสำคัญทางสาธารณสุขIntroduction and Public-Health Relevance

Local Moran's I (LISA) ตรวจสอบ cluster ระดับท้องที่ — แต่ละพื้นที่ i ได้ค่า I_i เฉพาะตัว, มี p-value จาก permutation test, และจัดประเภทเข้า 5 quadrants: HH (สูง-รายล้อมด้วยสูง = "hot cluster"), LL (ต่ำ-รายล้อมด้วยต่ำ = "cold cluster"), HL (สูงท่ามกลางต่ำ = "spatial outlier"), LH (ต่ำท่ามกลางสูง = "outlier"), และ NS (ไม่มีนัยสำคัญ). พัฒนาโดย Anselin (1995) เพื่อแก้ข้อจำกัดของ Moran's I ระดับ global ที่บอกได้แค่ว่ามี autocorrelation หรือไม่ แต่ไม่บอกว่าอยู่ที่ไหน.

Local Moran's I (LISA) detects local clusters — each location i gets its own I_i, a permutation p-value, and a 5-quadrant assignment: HH (high surrounded by high = "hot cluster"), LL (low-low cold cluster), HL (high among lows = "spatial outlier"), LH (low among highs = outlier), and NS (not significant). Anselin (1995) introduced LISA precisely to remedy the limitation of the global Moran's I, which only signals whether autocorrelation exists, not where.

2หลักการและสูตรLogic and Formula

กำหนดให้ z_i = ค่าของ x ที่ได้ standardise (หักด้วยค่าเฉลี่ย หารด้วย σ), และ w_ij = น้ำหนักเชิงพื้นที่ระหว่าง i กับ j (row-standardised):

Let z_i be the standardised value of x (centred & scaled by σ) and w_ij the row-standardised spatial weight between i and j:

p-value ของ I_i คำนวณจาก conditional permutation: สำหรับแต่ละ i, สุ่ม permute ค่า z_j ของพื้นที่อื่น 999 ครั้ง (default), คำนวณ I_i แต่ละครั้ง, แล้ววัดสัดส่วนที่ |I_i_permute| ≥ |I_i_observed|. การ permute เป็นแบบ conditional เพราะคงค่า z_i เดิมไว้.

p-value for I_i uses conditional permutation: for each i, permute the other z_j values 999 times (default), recompute I_i each time, and report the proportion of |I_i_permute| ≥ |I_i_observed|. The permutation is conditional because z_i is held fixed.

3สคริปต์ R และการทวนสอบกับ DAPHR Script and Verification against DAPH

library(spdep)
library(spData)

data(columbus)
col.gal.nb <- read.gal(system.file("etc/weights/columbus.gal", package = "spData"))
listw <- nb2listw(col.gal.nb, style = "W")

# Global Moran's I (for reference)
moran.test(columbus$CRIME, listw)
# I = 0.5117, p = 1.06e-09

# Local Moran's I
set.seed(12345)
loc <- localmoran(columbus$CRIME, listw, nsim = 999)
head(loc)        # → columns: Ii, E.Ii, Var.Ii, Z.Ii, Pr(z != E(Ii)) sim
# Cluster classification (HH/LL/HL/LH) via mean-centring + lag sign
ci   <- columbus$CRIME - mean(columbus$CRIME)
lagi <- lag.listw(listw, columbus$CRIME) - mean(columbus$CRIME)
quad <- ifelse(ci > 0 & lagi > 0, "HH",
        ifelse(ci < 0 & lagi < 0, "LL",
        ifelse(ci > 0 & lagi < 0, "HL", "LH")))
quad[loc[, "Pr(z != E(Ii)) sim"] > 0.05] <- "NS"
table(quad)

# Global Moran's I:
#   I = 0.51168     p-value = 1.06e-09 (analytical)

# Local Moran's I — top 5 by |Ii|:
#   POLYID  Ii        Z.Ii     p
#   1       1.16      2.78     0.005
#   3       0.59      1.71     0.045
#   13      0.84      2.31     0.018
#   17     -0.78     -2.04     0.041
#   34      0.96      2.56     0.011

# Quadrant assignment (Anselin 1995 convention):
#   HH = 6   LL = 11   HL = 1   LH = 2   NS = 29   (n = 49)

κ = 1.0 หมายถึง cluster assignment ตรงทุก location · "ความต่างสูงสุด" คือ ค่าผิดพลาดสูงสุดของ local I ที่ไหนสักแห่งใน dataset.κ = 1.0 means cluster labels agree at every location · "Max difference" is the maximum location-wise error of local I anywhere in the dataset.

4การแปลผลและตัวอย่างInterpretation and Example

4 quadrant ตาม Anselin 1995:
① HH (สีแดงเข้ม) — พื้นที่นี้และเพื่อนบ้านมีค่าสูงทั้งคู่ → hot cluster ที่ควรเฝ้าระวัง
② LL (สีน้ำเงินเข้ม) — ทั้งคู่ต่ำ → cold cluster
③ HL (สีชมพู) — ตัวเองสูง เพื่อนบ้านต่ำ → outlier ที่อาจเป็น "เกาะร้อนกลางพื้นที่ปลอดภัย" (ควรสอบสวน)
④ LH (สีฟ้าอ่อน) — ตัวเองต่ำ เพื่อนบ้านสูง → outlier (ทำไมที่นี่ปลอดภัย ทั้งที่รอบข้างไม่ปลอดภัย?)
⑤ NS (เทา) — ไม่มีนัยสำคัญที่ p < .05

The 4 Anselin (1995) quadrants:
① HH (dark red) — both this place and its neighbours are high → hot cluster worth surveillance
② LL (dark blue) — both low → cold cluster
③ HL (pink) — high surrounded by lows → outlier, possibly a "hot island in a safe area" (investigate)
④ LH (light blue) — low surrounded by highs → outlier (why is this one safe when surroundings are not?)
⑤ NS (grey) — not significant at p < .05

1บทนำและความสำคัญทางสาธารณสุขIntroduction and Public-Health Relevance

Bivariate Local Moran's I ขยาย LISA ให้ตรวจ ความสัมพันธ์เชิงพื้นที่ระหว่างตัวแปร 2 ตัว — ค่า x ที่พื้นที่ i กับค่า y ของเพื่อนบ้าน. ใช้แทน "correlation" แบบที่นับ across-unit ตรงไปตรงมาในกรณีที่ทั้งสองตัวแปรอาจมี spatial autocorrelation ของตัวเอง. พัฒนาโดย Anselin, Syabri & Smirnov (2002).

Bivariate Local Moran's I extends LISA to detect spatial co-variation between two variables — value of x at location i versus value of y in its neighbourhood. It replaces naïve cross-unit correlation when both variables have their own spatial autocorrelation. Introduced by Anselin, Syabri & Smirnov (2002).

2หลักการและสูตรLogic and Formula

กำหนดให้ z_i^x และ z_i^y เป็นค่า standardise ของ x และ y, w_ij = น้ำหนักเชิงพื้นที่:

Let z_i^x and z_i^y be standardised x and y, with w_ij the spatial weight:

p-value ใช้ permutation test เหมือน Univariate LISA — สุ่ม z_j^y 999 ครั้ง (คงค่า z_i^x ไว้) แล้ววัดสัดส่วน |I_perm| ≥ |I_obs|.

The p-value uses the same conditional permutation as Univariate LISA — permute z_j^y 999 times (holding z_i^x fixed) and measure the proportion of |I_perm| ≥ |I_obs|.

3สคริปต์ R และการทวนสอบกับ DAPHR Script and Verification against DAPH

หมายเหตุ: R package spdep ไม่มี Bivariate LISA สำเร็จรูป — ผู้ใช้ R สามารถคำนวณด้วยตนเองโดยใช้ spdep::lag.listw เพื่อหา spatial lag ของ y แล้วคูณกับ x ในรูปแบบ standardised values.

Note: R's spdep doesn't ship Bivariate LISA out of the box — R users can compute it manually with spdep::lag.listw for the spatial lag of y multiplied by standardised x.

library(spdep); library(spData)
data(columbus); listw <- nb2listw(read.gal(
  system.file("etc/weights/columbus.gal", package="spData")), style="W")

# Standardise CRIME and HOVAL
zx <- scale(columbus$CRIME)
zy <- scale(columbus$HOVAL)

# Lag of zy
lag_zy <- lag.listw(listw, zy)

# Bivariate Local Moran's I — observed
Ixy_obs <- zx * lag_zy

# Permutation p-value (1 location at a time, here we show vectorised)
set.seed(12345)
B <- 999
perm <- replicate(B, {
  zy_p <- sample(zy)            # permute y
  zx * lag.listw(listw, zy_p)   # bivariate I under H0
})
p <- rowMeans(abs(perm) >= abs(Ixy_obs))
head(data.frame(POLYID=1:length(zx), Ixy=Ixy_obs, p=p))

* ความแตกต่างเล็กน้อยจาก R manual เกิดจาก permutation seed ของ sample() ใน R ที่สุ่มต่างจาก seed=12345 ที่ DAPH ใช้ — ค่า I_i เองตรงทุกหลัก เพียงแต่จำนวน HH/LL cluster ที่ถูก flag significant อาจต่างกัน 1-2 จุดจากความ random ของ permutation.* The small gap from R manual is because R's sample() uses a different RNG path than DAPH's seed=12345 — the I_i values themselves match exactly; only 1–2 borderline-p locations may flip significance flag due to permutation randomness.

4การแปลผลและตัวอย่างInterpretation and Example

ใช้ quadrant แบบเดียวกับ Univariate แต่ความหมายต่างออกไป:
① HH: x สูงที่นี่ + y สูงรอบๆ → "x ที่นี่ขับเคลื่อน y ของเพื่อนบ้าน"
② LL: x ต่ำที่นี่ + y ต่ำรอบๆ
③ HL: x สูงที่นี่ + y ต่ำรอบๆ → ความขัดแย้ง (อาจเป็น factor ที่ไม่ได้สื่อให้เพื่อนบ้าน)
④ LH: x ต่ำที่นี่ + y สูงรอบๆ

Quadrants are the same as Univariate but mean different things:
① HH: high x here + high y around → "x here drives y in neighbours"
② LL: low x here + low y around
③ HL: high x here + low y around → discordant (maybe x doesn't reach neighbours)
④ LH: low x here + high y around

1บทนำและความสำคัญทางสาธารณสุขIntroduction and Public-Health Relevance

Getis-Ord statistics เป็นการตรวจ hotspot ตามความเข้มข้นของค่า (concentrations of value) — แตกต่างจาก Moran's I ที่ตรวจ ความคล้ายกัน ทั้งสูงและต่ำ. มี 3 สถิติคู่กัน: Global G (รวมทั้ง dataset), Local G_i (เฉพาะ neighbour ของ i), และ Local G*_i (รวม i เองด้วย). พัฒนาโดย Getis & Ord (1992) และ Ord & Getis (1995). ใช้แพร่หลายในการระบุจุดร้อน/จุดเย็นของโรคใน epidemiology — ค่า G ที่ใหญ่และเป็นบวก = hotspot, ใหญ่และเป็นลบ = coldspot.

Getis-Ord statistics detect concentrations of value magnitude — distinct from Moran's I which detects similarity (both high or both low). The trio: Global G (whole dataset), Local G_i (neighbours of i), and Local G*_i (neighbours plus i). Introduced by Getis & Ord (1992) and Ord & Getis (1995). Widely used in disease epidemiology to identify hotspots/coldspots — large positive G = hotspot, large negative G = coldspot.

2หลักการและสูตรLogic and Formula

ค่า G แปลความหมายผ่าน z-score:

Significance is read via the z-score:

3สคริปต์ R และการทวนสอบกับ DAPHR Script and Verification against DAPH

library(spdep); library(spData)
data(columbus)
nb <- read.gal(system.file("etc/weights/columbus.gal", package="spData"))

# Distance-based weights for Getis-Ord (queen contiguity also works)
listw  <- nb2listw(nb, style = "B")        # binary
listws <- nb2listw(include.self(nb), style="B")  # for G* (self-included)

# Global G
globalG.test(columbus$CRIME, listw)
# expect: G = 0.0098, z = 1.86, p = 0.031 (one-sided)

# Local G_i  (excludes i)
gi  <- localG(columbus$CRIME, listw)
# Local G*_i (includes i)
gis <- localG(columbus$CRIME, listws)
head(data.frame(POLYID=1:length(gi), Gi=as.numeric(gi), Gi_star=as.numeric(gis)))

# Global G:    G = 0.00978,   z = 1.864,    p = 0.031 (one-sided positive)
# Local G_i — top hotspots (z > 1.96):
#   POLYID 1    Gi*=2.31,   z=2.74,  p=0.006  → hotspot
#   POLYID 34   Gi*=1.85,   z=2.31,  p=0.010  → hotspot
#   POLYID 17   Gi*=-1.74,  z=-2.05, p=0.040  → coldspot

4การแปลผลและตัวอย่างInterpretation and Example

การเลือกใช้ G_i เทียบกับ G*_i:
— G_i ตรวจว่า "เพื่อนบ้านของ i มีค่าสูงผิดปกติไหม" — ตัวเองไม่นับ
— G*_i ตรวจว่า "พื้นที่ i รวมเพื่อนบ้านมีค่าสูงผิดปกติไหม" — รวมตัวเอง
โดยทั่วไป G*_i ใช้บ่อยกว่า เพราะตอบคำถาม policy ตรงกว่า ("ตำบลนี้เป็น hotspot ไหม") · DAPH ออกผลทั้งคู่ให้.

G_i vs G*_i — which to use:
— G_i: "are i's neighbours unusually high?" — excludes i itself
— G*_i: "is the i+neighbours area unusually high?" — includes i
G*_i is preferred in most policy applications ("is this sub-district a hotspot?"). DAPH reports both.

1ทำไมการเลือก weights ถึงสำคัญที่สุดWhy the weights choice matters most

ทุกสถิติในบทนี้เริ่มจากเมทริกซ์ W ที่บอกว่า "พื้นที่ไหนเป็นเพื่อนบ้านกับพื้นที่ไหน" · ถ้า W เปลี่ยน ค่า ρ, λ, Moran's I และข้อสรุปก็เปลี่ยนตาม ดังนั้น ต้องรายงานชนิดของ weights ที่ใช้เสมอ — DAPH พิมพ์บรรทัดนี้ให้อัตโนมัติทั้งใน Tab3 และไฟล์ Word

Every statistic in this chapter starts from a matrix W that declares which areas are neighbours. Change W and ρ, λ, Moran's I and the conclusion change with it, so the weights specification must always be reported — DAPH prints that line automatically in Tab3 and in the Word file.

2ตัวเลือกที่มีให้ และควรใช้เมื่อไรThe options, and when to use each

3คำเตือน "เชื่อมกันไม่ครบ" หมายความว่าอย่างไรWhat the "not fully connected" warning means

ถ้า weights ที่เลือกทำให้แผนที่แตกออกเป็นหลายก้อนที่ไม่มีเพื่อนบ้านร่วมกัน DAPH จะต่อท้ายบรรทัด spatial weights ว่า "⚠ เชื่อมกันไม่ครบ: แยกเป็น N กลุ่ม…" · กรณีนี้เมทริกซ์ W จะมีค่า eigenvalue เท่ากับ 1 หลายค่า ทำให้ พารามิเตอร์เชิงพื้นที่ถูกระบุค่าได้อ่อน — ค่า ρ ที่ได้อาจต่างจากโปรแกรมอื่นเกิน 10% และตัวประมาณค่าแบบ RE-Error อาจไม่ลู่เข้าเลย · ทางแก้: เปลี่ยนเป็น "Queen ปรับ" (ถ้าสาเหตุคือเกาะ) หรือขยายรัศมี หรือใช้ KNN

If the chosen weights split the map into groups with no shared neighbour, DAPH appends "⚠ not fully connected: N separate groups…" to the weights line. W then carries a unit eigenvalue per group, so the spatial parameter is weakly identified — ρ can differ from another program by more than 10%, and the RE-Error estimator may fail to converge. Fix: switch to the modified Queen (if islands are the cause), widen the radius, or use KNN.

ถ้าเกาะอยู่ไกลเกินกว่าจะถือเป็นเพื่อนบ้านจริง ให้กรอกรัศมีสูงสุดที่ยอมรับได้ในช่อง threshold — เกาะที่ไกลกว่านั้นจะไม่ถูกเชื่อม และโปรแกรมจะระบุชื่อพร้อมระยะทางไว้ในรายงาน

If an island is genuinely too far to count as a neighbour, enter the largest acceptable radius in the threshold box — islands beyond it stay unlinked and the report names them with their distance.

4เอกสารอ้างอิงReferences

• Anselin, L. (1988). Spatial Econometrics: Methods and Models. Dordrecht: Kluwer Academic Publishers.
• Getis, A., & Aldstadt, J. (2004). Constructing the spatial weights matrix using a local statistic. Geographical Analysis, 36(2), 90–104. doi:10.1111/j.1538-4632.2004.tb01127.x

1บทนำและความสำคัญทางสาธารณสุขIntroduction and Public-Health Relevance

Distribution Map ใน DAPH รวม 2 เทคนิคการ visualize ตัวแปรเชิงพื้นที่: Choropleth (แบ่งสีตามคลาส) และ KDE Heatmap (Kernel Density Estimation บนพิกัด point data). Choropleth ใช้แบ่ง polygon เป็นกลุ่มสี — DAPH รองรับ 4 วิธีหลัก: Quantile, Equal Interval, Natural Breaks (Jenks 1967), และ Standard Deviation. KDE Heatmap ใช้ smooth point data ออกเป็น density surface ตามสูตร Silverman (1986).

Distribution Map in DAPH bundles two spatial-visualization techniques: Choropleth (colour by class) and KDE Heatmap (Kernel Density Estimation on point data). Choropleth partitions polygons into colour classes — DAPH supports 4 schemes: Quantile, Equal Interval, Natural Breaks (Jenks 1967), and Standard Deviation. KDE Heatmap smooths point data into a density surface per Silverman (1986).

2หลักการและสูตรLogic and Formula

3สคริปต์ R และการทวนสอบกับ DAPHR Script and Verification against DAPH

library(classInt); library(spData); library(KernSmooth)
data(columbus)
x <- columbus$CRIME

# Quantile (5 classes)
ci_q  <- classIntervals(x, n = 5, style = "quantile")
# Equal Interval
ci_e  <- classIntervals(x, n = 5, style = "equal")
# Natural Breaks (Jenks 1967)
ci_j  <- classIntervals(x, n = 5, style = "jenks")
# Standard Deviation
ci_sd <- classIntervals(x, n = 5, style = "sd")

print(ci_q$brks);  print(ci_e$brks);  print(ci_j$brks);  print(ci_sd$brks)

# KDE on point data (e.g. CRIME centroid coordinates)
xy <- coordinates(columbus[, c("X","Y")])
kde <- bkde2D(xy, bandwidth = c(dpik(xy[,1]), dpik(xy[,2])))

# Quantile:        0.18   15.4   26.0   38.2   54.8   68.9
# Equal Interval:  0.18   13.9   27.7   41.4   55.2   68.9
# Jenks:           0.18   15.0   28.4   42.4   55.6   68.9
# Std Dev (μ=35.1, σ=16.7):
#                  0.18   18.4   26.7   35.1   43.4   51.7   68.9   (7 breaks for ±2σ scheme)

Jenks algorithm ของ DAPH ใช้สูตร Fisher-Jenks (1958) ที่เป็นรูป optimal ทางคณิตศาสตร์ — ผลตรงกับ R classInt v0.4-x ทุกหลัก.DAPH's Jenks uses the Fisher-Jenks (1958) algorithm which is mathematically optimal — exact agreement with R classInt v0.4-x.

4การแปลผลและตัวอย่างInterpretation and Example

เลือกวิธี classify อย่างไร:
— Quantile — แต่ละ class มีจำนวน area เท่ากัน · เหมาะกับการเปรียบเทียบลำดับ (rank)
— Equal Interval — ช่วงค่าตัดเท่ากัน · เหมาะเมื่อข้อมูลกระจายเกือบ uniform
— Natural Breaks (Jenks) — ลด within-class variance · เหมาะกับ data ที่มี cluster ตามธรรมชาติ (default ของ DAPH)
— Standard Deviation — แบ่งตาม ±σ จากค่าเฉลี่ย · เหมาะเมื่อข้อมูลเป็น Normal และต้องการเน้น outlier

Choosing a classification scheme:
— Quantile — equal count per class · best for ranking comparisons
— Equal Interval — equal range per class · best when data are roughly uniform
— Natural Breaks (Jenks) — minimises within-class variance · best when data have natural clusters (DAPH default)
— Standard Deviation — split by ±σ around the mean · best for Normally distributed data emphasising outliers