预测核辐射强度的分布

空间数据分析

黄湘云

统计之都 @ 北京会议中心

2023年7月12日

目录

  1. 核污染研究的背景和意义

  2. 数据探索和分析的过程

    • 描述数据的空间分布
    • 广义线性模型及模型诊断
    • 空间广义线性混合效应模型
    • 模型拟合及准备预测位置
    • 模型预测及展示预测分布

氢弹核试验

这颗氢弹的当量相当于广岛原子弹的 1000 倍。 由于核辐射扩散,导致朗格拉普环礁的主岛受到核辐射1

朗格拉普岛

图 1: 朗格拉普岛地理位置特殊,是一个重要的军事基地,建有环岛公路和机场。卫星拍摄的影像数据(2020 年)来自 https://rmi-data.sprep.org/

采集的核辐射数据(部分)

表格 1: 采集的核辐射数据

(a) 核辐射检测数据
横坐标 纵坐标 粒子数目 统计时间
-6050 -3270 75 300
-6050 -3165 371 300
-5925 -3320 1931 300
-5925 -3165 4357 300
-5800 -3350 2114 300
-5800 -3165 2318 300
(b) 海岸线坐标数据
横坐标 纵坐标
-5509.236 -3577.438
-5544.821 -3582.250
-5561.604 -3576.926
-5580.780 -3574.535
-5599.687 -3564.288
-5605.922 -3560.910

采样数据的空间分布

图 2: 采样点在岛上的分布

图 3: 采样点在岛上的分布

岛上各采样点的核辐射强度

图 4: 岛上各采样点的核辐射强度

广义线性模型及拟合

核辐射是由放射元素衰变产生的,通常用单位时间释放出来的粒子数目表示辐射强度,因此,建立如下泊松型广义线性模型来拟合核辐射强度。

\begin{aligned} \log(\lambda_i) &= \beta \\ y_i & \sim \mathrm{Poisson}(t_i\lambda_i) \end{aligned}

待估参数 \beta

fit_rongelap_poisson <- glm(counts ~ 1,
  family = poisson(link = "log"), 
  offset = log(time), data = rongelap
)
summary(fit_rongelap_poisson)

Call:
glm(formula = counts ~ 1, family = poisson(link = "log"), data = rongelap, 
    offset = log(time))

Coefficients:
            Estimate Std. Error z value Pr(>|z|)    
(Intercept) 2.013954   0.001454    1385   <2e-16 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

(Dispersion parameter for poisson family taken to be 1)

    Null deviance: 61567  on 156  degrees of freedom
Residual deviance: 61567  on 156  degrees of freedom
AIC: 63089

Number of Fisher Scoring iterations: 4

模型诊断 I

图 5: 残差的空间分布

颜色深和颜色浅的点分别聚集在一起,且与周围点的颜色呈现层次变化,拟合残差存在一定的空间相关性。

模型诊断 II