fix bug when pred variance = 0

swpease · Sep 15, 2023 · 32a6d5e · 32a6d5e
1 parent 941a78b
commit 32a6d5e
Show file tree

Hide file tree

Showing 3 changed files with 6 additions and 1 deletion.
diff --git a/R/compute_edf.R b/R/compute_edf.R
@@ -73,11 +73,16 @@ compute_edf = function(mgcv_model, object, rho_names, sigma_raw_names){
 
     # Because some mvgam families do not (yet) have simple functions to compute
     # prediction variance, we will resort to brute force sampling
-    preds <- do.call(rbind, lapply(1:100, function(x){
+    preds <- do.call(rbind, lapply(1:75, function(x){
       # luckily the predict function is vectorized!!
       predict(object, process_error = FALSE, type = 'response')[best_draw,]
     }))
     pred_variance <- apply(preds, 2, var)
+    if(any(pred_variance == 0)){
+      pred_variance[which(pred_variance == 0)] <-
+        mu[which(pred_variance == 0)]
+    }
+
     w <- as.numeric(mgcv_model$family$mu.eta(eta) / pred_variance)
     XWX <- t(X) %*% (w * X)
     rho <- mgcv_model$sp

diff --git a/src/mvgam.dll b/src/mvgam.dll
diff --git a/tests/testthat/Rplots.pdf b/tests/testthat/Rplots.pdf