使用R语言moveHMM包绘制HMM模型时出现逻辑值缺失错误

问题背景 #

使用moveHMM包拟合3状态隐马尔可夫模型，代码如下：

hmm_model <- fitHMM(data = move_data, nbStates = 3, stepPar0 = stepPar0, anglePar0 = anglePar0, formula = ~1)

模型结果符合生物学意义且无警告，但执行plot(hmm_model, type = "states")或plot(hmm_model)时触发错误：

> plot(hmm_model, type = "states")
Decoding states sequence... DONE
Error in if (maxdens > ymax & maxdens < 1.5 * ymax) { : 
  missing value where TRUE/FALSE needed

模型结构信息 #

> str(hmm_model)
List of 7
 $ data       :'data.frame':    16515 obs. of  6 variables:
  ..$ ID         : chr [1:16515] "1" "1" "1" "1" ...
  ..$ step       : num [1:16515] 0.327 0.346 0.346 0.346 0.346 ...
  ..$ angle      : num [1:16515] NA -0.156381 -0.004527 -0.001887 -0.000926 ...
  ..$ x          : num [1:16515] 582629 582629 582628 582628 582628 ...
  ..$ y          : num [1:16515] -2195359 -2195359 -2195358 -2195358 -2195358 ...
  ..$ (Intercept): num [1:16515] 1 1 1 1 1 1 1 1 1 1 ...
 $ mle        :List of 5
  ..$ stepPar : num [1:3, 1:3] 3.43e-15 8.36e+09 1.00 8.31e-02 1.13e-01 ...
  .. ..- attr(*, "dimnames")=List of 2
  .. .. ..$ : chr [1:3] "mean" "sd" "zero-mass"
  .. .. ..$ : chr [1:3] "state 1" "state 2" "state 3"
  ..$ anglePar: num [1:2, 1:3] 1.31e-03 2.83 1.32e-03 2.23e+01 -3.31e-10 ...
  .. ..- attr(*, "dimnames")=List of 2
  .. .. ..$ : chr [1:2] "mean" "concentration"
  .. .. ..$ : chr [1:3] "state 1" "state 2" "state 3"
  ..$ beta    : num [1, 1:6] -6.91 -37.97 -9.65 -25.74 -0.43 ...
  .. ..- attr(*, "dimnames")=List of 2
  .. .. ..$ : chr "intercept"
  .. .. ..$ : chr [1:6] "1 -> 2" "1 -> 3" "2 -> 1" "2 -> 3" ...
  ..$ delta   : num [1:3] 2.67e-10 1.00 1.32e-14
  ..$ gamma   : num [1:3, 1:3] 9.99e-01 6.44e-05 3.38e-01 9.93e-04 1.00 ...
 $ mod        :List of 6
  ..$ minimum   : num -27526
  ..$ estimate  : num [1:23] -33.31 -2.49 1.61 22.85 -2.18 ...
  ..$ gradient  : num [1:23] 0 -0.00268 0 0 0.00216 ...
  ..$ hessian   : num [1:23, 1:23] 0 0 0 0 0 0 0 0 0 0 ...
  ..$ code      : int 1
  ..$ iterations: int 73
 $ conditions :List of 6
  ..$ stepDist     : chr "gamma"
  ..$ angleDist    : chr "vm"
  ..$ zeroInflation: logi TRUE
  ..$ estAngleMean : logi TRUE
  ..$ stationary   : logi FALSE
  ..$ formula      :Class 'formula'  language ~1
  .. .. ..- attr(*, ".Environment")=<environment: R_GlobalEnv> 
 $ rawCovs    :'data.frame':    16515 obs. of  0 variables
 $ knownStates: NULL
 $ nlmTime    : 'proc_time' Named num [1:5] 66.16 0.72 66.88 NA NA
  ..- attr(*, "names")= chr [1:5] "user.self" "sys.self" "elapsed" "user.child" ...
 - attr(*, "class")= chr [1:2] "moveHMM" "list"

解决建议 #

• 检查极端参数值：从模型参数看，state2的step均值为8.36e+09，远大于数据中step的实际值（数据中step最大约0.346），这会导致计算该状态的step密度时出现无穷大或NA，触发绘图判断条件的错误。先确认数据中是否存在异常大的step值，若没有则说明拟合陷入局部最优。
• 手动验证密度计算：用moveHMM内置函数验证各状态的密度是否正常，示例代码：

  # 测试state2的step密度
  x_seq <- seq(0, max(move_data$step, na.rm=TRUE), length.out=100)
  state2_step_dens <- dstep(x_seq, dist="gamma", par=hmm_model$mle$stepPar[,2], zeroInflation=TRUE)
  # 检查是否有NA或无穷大
  any(is.na(state2_step_dens))
  any(is.infinite(state2_step_dens))
  

  如果返回TRUE，说明该状态的密度计算异常，是报错的直接原因。
• 重新设置初始值拟合：当前异常参数大概率是初始值stepPar0设置不合理导致的。建议根据数据分位数设置初始值：比如将step数据按大小分成3组，取每组的均值、标准差和零比例作为对应状态的初始参数，再重新拟合模型。
• 排查零膨胀设置：尝试拟合不带零膨胀的模型（zeroInflation=FALSE），若绘图正常，则说明零膨胀参数的拟合存在问题，需调整初始值或检查数据中零步的比例分布。
• 更新moveHMM版本：部分旧版本可能存在绘图函数的bug，更新到最新版包后重试。
• 检查数据缺失值：虽然moveHMM会处理angle列的缺失值，但确认数据中是否存在其他异常值（如step为负、极端大值等），清理后再重新拟合。

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