如何在Hmsc中正确指定日期嵌套的结构化时间随机效应？

在两块森林样地开展7天熊蜂群落监测，监测时段为日出至日落，每日设置16个时间点（日出/日落仅在其中一块样地监测，其余时间点两块样地均采样）。需在Hmsc模型中定义两类随机效应：

• 样地水平的随机效应
• 嵌套在日期内的结构化时间随机效应，用于模拟单日监测中逐小时观测的群落相关性

现有数据结构及建模代码如下：

数据结构 #

head(XData)
date  plot time flowers    av.li    av.tp duration
1 220512  open    0     8.4 4.615006 6.861538       60
2 220512  open    1     8.4 6.973745 7.333333       15
3 220512 dense    1     6.6 7.387100 7.666667       15
4 220512  open    2     8.4 7.819142 7.975000       15
5 220512 dense    2     6.6 8.082691 8.040000       15
6 220512  open    3     8.4 8.854704 9.040000       15

现有建模代码 #

研究设计矩阵定义 #

# Definition of the study design matrix -
studyDesign <- matrix(NA, nrow(data), 3)
studyDesign <- as.data.frame(studyDesign)
studyDesign[,1] <- XData$plot
studyDesign[,1] <- as.factor(studyDesign[,1])
studyDesign[,2] <- paste('time', XData$time, sep = "_")
studyDesign[,2] <- as.factor(studyDesign[,2])
studyDesign[,3] <- paste('date', XData$date, sep = "_")
studyDesign[,3] <- as.factor(studyDesign[,3])
colnames(studyDesign) <- c("plot", "time", "date")

head(studyDesign)
   plot   time        date
1  open time_0 date_220512
2  open time_1 date_220512
3 dense time_1 date_220512
4  open time_2 date_220512
5 dense time_2 date_220512
6  open time_3 date_220512

注：也曾将time和date分别定义为0（日出）到17（日落）的因子，date格式为220512（代表2022年5月12日）。

随机效应与模型定义 #

# TIME
time <- unique(studyDesign[,2])
xy <- as.data.frame(XData %>% dplyr::select(time))
xy <- unique(xy)
colnames(xy) <- c("time")
sRL <- xy
rownames(sRL) <- time
rL_T <- HmscRandomLevel(sData = sRL)

# PLOT AND DATE
rL2 <- HmscRandomLevel(units = unique(studyDesign$plot))
rL3 <- HmscRandomLevel(units = unique(studyDesign$date))

# DEFINITION OF THE MODEL
m.r <- Hmsc(Y = Y1,
            XData = XData, XFormula = XFormula,
            TrData = TrData, TrFormula = TrFormula,
            distr = "probit",
            ranLevels = list("plot"=rL2, "date"=rL3,"time"=rL_T ),
            studyDesign = studyDesign)

模型可正常运行，但无法确认是否正确指定了重复采样的时间结构，尤其担心时间因子被错误排序（例如按1,11,12…17,2,3…9排序而非0,1,2…17），手动排序因子后结果仍未改善，需指导如何正确指定嵌套的结构化时间随机效应并保证时间顺序正确。

1. 修复时间因子排序问题 #

当前代码中用paste('time', XData$time, sep = "_")生成的时间因子会按字符排序（如time_10排在time_2前面），需手动指定因子水平的正确顺序：

# 基于原始time数值生成有序因子，确保时间点按0→17排列
studyDesign$time <- factor(paste('time', XData$time, sep = "_"),
                          levels = paste('time', 0:17, sep = "_"),
                          ordered = TRUE)

# 若直接使用原始time数值作为因子，同样指定水平顺序
# studyDesign$time <- factor(XData$time, levels = 0:17, ordered = TRUE)

2. 正确指定嵌套式结构化时间随机效应 #

原设置中time是全局随机水平，但实际需求是每个日期内的时间点为独立嵌套结构（不同日期的同一时间点无需共享结构），需按以下步骤调整：

步骤1：构建嵌套研究设计矩阵 #

创建date_time复合标识，代表每个日期下的独立时间单元：

# 生成日期-时间复合因子，明确嵌套关系
studyDesign$date_time <- paste(studyDesign$date, studyDesign$time, sep = "_")
studyDesign$date_time <- factor(studyDesign$date_time)

步骤2：为嵌套时间单元构建结构化协变量 #

Hmsc的结构化随机效应依赖连续时间协变量建模相关性，需为每个date_time提取对应时间数值：

# 提取每个date_time对应的时间数值（0-17）
time_vals <- as.numeric(gsub("time_", "", studyDesign$time))
# 按date_time分组，保留唯一的时间协变量数据
sRL_data <- data.frame(time_val = time_vals)
rownames(sRL_data) <- studyDesign$date_time
sRL_data <- unique(sRL_data)

# 创建结构化随机水平，基于时间数值构建相关性矩阵
rL_date_time <- HmscRandomLevel(sData = sRL_data)

步骤3：重新定义模型随机效应 #

将嵌套的date_time作为结构化随机水平，保留plot作为非结构化随机水平：

m.r <- Hmsc(Y = Y1,
            XData = XData, XFormula = XFormula,
            TrData = TrData, TrFormula = TrFormula,
            distr = "probit",
            ranLevels = list("plot" = rL2, "date_time" = rL_date_time),
            studyDesign = studyDesign)

3. 验证时间结构正确性 #

模型运行后，可提取结构化随机效应的相关性矩阵验证：

# 提取后验估计中的随机效应相关性矩阵
posterior <- getPostEstimate(m.r)
corr_matrix <- posterior$RanLevels$date_time$Cor[[1]]
# 查看前5个时间单元的相关性，相邻时间点相关性应更高
print(corr_matrix[1:5, 1:5])

内容的提问来源于stack exchange，提问作者OceBartho