From e073bf3ff2572bd472a7bb3d7afbfed37a537588 Mon Sep 17 00:00:00 2001 From: Laura Marshall Date: Mon, 27 Jan 2025 13:23:48 +0000 Subject: [PATCH 01/14] Issue #96 gof dfs for MR models fixes when there are insufficient df for test for io models --- DESCRIPTION | 2 +- NEWS.md | 1 + R/gof.io.R | 17 ++++-- R/print.ddf.gof.R | 17 ++++-- tests/testthat/test_gof.R | 117 ++++++++++++++++++++++++++++++++++++++ 5 files changed, 143 insertions(+), 11 deletions(-) create mode 100644 tests/testthat/test_gof.R diff --git a/DESCRIPTION b/DESCRIPTION index 9709ad72..0d15b210 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -18,7 +18,7 @@ Description: Animal abundance estimation via conventional, multiple covariate fitting is performed via maximum likelihood. Also included are diagnostics and plotting for fitted detection functions. Abundance estimation is via a Horvitz-Thompson-like estimator. -Version: 3.0.0.9001 +Version: 3.0.0.9002 URL: https://github.com/DistanceDevelopment/mrds/ BugReports: https://github.com/DistanceDevelopment/mrds/issues Depends: diff --git a/NEWS.md b/NEWS.md index 231f16b3..d18495a8 100644 --- a/NEWS.md +++ b/NEWS.md @@ -3,6 +3,7 @@ Bug Fixes * Fixed formatting issue in flnl.grad help +* Fixed output for gof test when degrees of freedom for MR io models is <= 0 # mrds 3.0.0 diff --git a/R/gof.io.R b/R/gof.io.R index ad4a2b34..168f7e01 100644 --- a/R/gof.io.R +++ b/R/gof.io.R @@ -84,7 +84,17 @@ gof.io <- function(model, breaks=NULL, nc=NULL){ }else{ p.2 <- 1-pchisq(chisq.2,df.2) } - + + # Calculate the pooled chi-square + df.pool <- 3*nc-length(model$par)-1 + + if(df.pool <= 0){ + df.pool <- NA + p.pool <- NA + }else{ + p.pool <- 1-pchisq(chisq.1+chisq.2, df.pool) + } + return(list(chi1=list(observed=observed.count.1, expected=expected.1, chisq=chisq.1, @@ -96,7 +106,6 @@ gof.io <- function(model, breaks=NULL, nc=NULL){ p=p.2, df=df.2), pooled.chi=list(chisq=chisq.1+chisq.2, - df=3*nc-length(model$par)-1, - p=1-pchisq(chisq.1+chisq.2, - 3*nc-length(model$par)-1)))) + df=df.pool, + p=p.pool))) } diff --git a/R/print.ddf.gof.R b/R/print.ddf.gof.R index 2bb5005b..4bbaaf7d 100644 --- a/R/print.ddf.gof.R +++ b/R/print.ddf.gof.R @@ -83,12 +83,17 @@ print.ddf.gof <- function(x, digits=3, ...){ " with ", gof$chisquare$chi2$df, " degrees of freedom\n", sep="")) } - - cat(paste("\n\nTotal chi-square = ", - format(gof$chisquare$pooled.chi$chisq, digits=5), - " P = ", format(gof$chisquare$pooled.chi$p, digits=5), - " with ", gof$chisquare$pooled.chi$df, " degrees of freedom\n", - sep="")) + + if(!is.na(gof$chisquare$pooled.chi$p)){ + cat(paste("\n\nTotal chi-square = ", + format(gof$chisquare$pooled.chi$chisq, digits=5), + " P = ", format(gof$chisquare$pooled.chi$p, digits=5), + " with ", gof$chisquare$pooled.chi$df, " degrees of freedom\n", + sep="")) + }else{ + cat("\nTotal chi-square: No degrees of freedom for test\n") + } + } } diff --git a/tests/testthat/test_gof.R b/tests/testthat/test_gof.R new file mode 100644 index 00000000..7bc3f59d --- /dev/null +++ b/tests/testthat/test_gof.R @@ -0,0 +1,117 @@ +# Create test data (this was simulated using dsims but for time and not having this as a dependency the data are pasted in here). + +tot.data <- data.frame(Region.Label = "region", + Area = 1e+06, + Sample.Label = c(17, 17, 19, 19, 10, 10, 19, 19, 17, 17, 6, 6, 15, 15, 6, 6, 19, 19, 9, 9, 6, 6, 6, 6, 2, 2, 14, 14, 8, 8, 6, 6, 6, 6, 6, 6, 2, 2, 12, 12, 12, 12, 16, 16, 11, 11, 2, 2, 3, 3, 3, 3, 16, 16, 11, 11, 16, 16, 20, 20, 12, 12, 14, 14, 8, 8, 4, 4, 16, 16, 13, 13, 10, 10, 2, 2, 15, 15, 14, 14, 9, 9, 15, 15, 19, 19, 19, 19, 18, 18, 1, 1, 7, 7, 2, 2, 7, 7, 14, 14, 7, 7, 4, 4, 13, 13, 15, 15, 16, 16, 17, 17, 2, 2, 9, 9, 8, 8, 12, 12, 13, 13, 5, 5, 6, 6, 1, 1, 5, 5, 10, 10, 6, 6, 9, 9, 20, 20, 20, 20, 8, 8, 11, 11, 1, 1, 15, 15, 5, 5, 4, 4, 16, 16, 13, 13, 14, 14, 5, 5, 13, 13, 9, 9, 11, 11, 16, 16, 8, 8, 13, 13, 10, 10, 7, 7, 1, 1, 9, 9, 6, 6, 14, 14, 1, 1, 7, 7, 3, 3, 17, 17, 19, 19, 13, 13, 5, 5, 15, 15, 10, 10, 14, 14, 16, 16, 2, 2, 2, 2, 18, 18, 17, 17, 17, 17, 10, 10, 5, 5, 1, 1, 18, 18, 13, 13, 11, 11, 9, 9, 9, 9, 3, 3, 6, 6, 12, 12, 18, 18, 19, 19, 5, 5, 11, 11, 2, 2, 8, 8, 17, 17, 1, 1, 6, 6, 11, 11, 11, 11, 14, 14, 9, 9, 14, 14, 9, 9, 2, 2, 17, 17, 15, 15, 7, 7, 8, 8, 6, 6, 10, 10, 20, 20, 9, 9, 14, 14, 2, 2, 15, 15, 16, 16, 16, 16, 8, 8, 6, 6, 11, 11, 16, 16, 10, 10, 9, 9, 16, 16, 4, 4, 11, 11, 18, 18, 3, 3, 14, 14, 7, 7, 15, 15, 5, 5, 15, 15, 2, 2, 3, 3, 10, 10, 1, 1, 9, 9, 16, 16, 8, 8, 8, 8, 19, 19, 13, 13, 18, 18, 18, 18, 15, 15, 4, 4, 11, 11, 15, 15, 8, 8, 3, 3, 7, 7, 20, 20, 1, 1, 13, 13), + Effort = 500, + object = c(3, 3, 4, 4, 6, 6, 7, 7, 9, 9, 10, 10, 11, 11, 12, 12, 13, 13, 14, 14, 15, 15, 16, 16, 17, 17, 18, 18, 19, 19, 20, 20, 21, 21, 22, 22, 23, 23, 25, 25, 26, 26, 27, 27, 28, 28, 29, 29, 30, 30, 31, 31, 33, 33, 34, 34, 35, 35, 36, 36, 40, 40, 41, 41, 43, 43, 46, 46, 47, 47, 48, 48, 50, 50, 53, 53, 54, 54, 55, 55, 56, 56, 57, 57, 58, 58, 59, 59, 60, 60, 63, 63, 64, 64, 65, 65, 66, 66, 67, 67, 68, 68, 69, 69, 70, 70, 71, 71, 72, 72, 74, 74, 75, 75, 76, 76, 77, 77, 78, 78, 81, 81, 82, 82, 83, 83, 84, 84, 85, 85, 86, 86, 87, 87, 89, 89, 90, 90, 92, 92, 93, 93, 94, 94, 96, 96, 97, 97, 98, 98, 99, 99, 100, 100, 102, 102, 103, 103, 104, 104, 106, 106, 107, 107, 110, 110, 111, 111, 112, 112, 113, 113, 114, 114, 115, 115, 117, 117, 118, 118, 119, 119, 120, 120, 121, 121, 122, 122, 124, 124, 125, 125, 128, 128, 129, 129, 131, 131, 132, 132, 133, 133, 134, 134, 135, 135, 136, 136, 137, 137, 139, 139, 140, 140, 141, 141, 142, 142, 144, 144, 145, 145, 146, 146, 148, 148, 149, 149, 150, 150, 154, 154, 155, 155, 156, 156, 157, 157, 158, 158, 159, 159, 160, 160, 161, 161, 162, 162, 163, 163, 164, 164, 165, 165, 166, 166, 168, 168, 169, 169, 170, 170, 172, 172, 173, 173, 175, 175, 176, 176, 177, 177, 178, 178, 179, 179, 180, 180, 182, 182, 184, 184, 185, 185, 187, 187, 188, 188, 190, 190, 192, 192, 193, 193, 196, 196, 197, 197, 198, 198, 201, 201, 202, 202, 203, 203, 204, 204, 205, 205, 206, 206, 207, 207, 208, 208, 209, 209, 210, 210, 212, 212, 213, 213, 216, 216, 217, 217, 218, 218, 219, 219, 220, 220, 222, 222, 223, 223, 224, 224, 225, 225, 226, 226, 227, 227, 229, 229, 230, 230, 231, 231, 232, 232, 234, 234, 236, 236, 240, 240, 241, 241, 242, 242, 245, 245, 246, 246, 247, 247, 249, 249), + observer = rep(c(1,2),186), + detected = c(1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 0, 1, 1, 1, 1, 1, 0, 1, 1, 1, 0, 1, 1, 1, 0, 1, 1, 1, 1, 0, 0, 1, 1, 1, 1, 0, 0, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 0, 1, 1, 0, 1, 1, 1, 1, 1, 0, 1, 1, 1, 0, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 0, 1, 1, 0, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 0, 1, 0, 1, 0, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 0, 1, 0, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 0, 1, 1, 1, 0, 1, 1, 1, 1, 0, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 0, 1, 1, 1, 1, 1, 0, 1, 0, 1, 1, 0, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 0, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 0, 1, 0, 1, 1, 1, 1, 0, 0, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 0, 0, 1, 1, 0, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1), + distance = c(3.35, 3.35, 39.69, 39.69, 8.52, 8.52, 2.57, 2.57, 1.64, 1.64, 22.18, 22.18, 18.27, 18.27, 15.54, 15.54, 9.96, 9.96, 7.16, 7.16, 12.35, 12.35, 11.87, 11.87, 37.13, 37.13, 9.06, 9.06, 4.41, 4.41, 20.85, 20.85, 44.76, 44.76, 30, 30, 25.58, 25.58, 8.96, 8.96, 18.39, 18.39, 0.63, 0.63, 12.05, 12.05, 11.84, 11.84, 11.17, 11.17, 40.31, 40.31, 15.69, 15.69, 35.18, 35.18, 14.37, 14.37, 38.38, 38.38, 17.19, 17.19, 30.3, 30.3, 30.41, 30.41, 3.74, 3.74, 27.23, 27.23, 18.79, 18.79, 14.21, 14.21, 6, 6, 19.58, 19.58, 13.28, 13.28, 24.97, 24.97, 2.64, 2.64, 9.89, 9.89, 11.38, 11.38, 28.69, 28.69, 21.62, 21.62, 41.23, 41.23, 11.99, 11.99, 6.27, 6.27, 34.5, 34.5, 49.55, 49.55, 7.33, 7.33, 17.52, 17.52, 36.47, 36.47, 15.05, 15.05, 20.37, 20.37, 7.08, 7.08, 29.22, 29.22, 14.57, 14.57, 32.39, 32.39, 27, 27, 48.06, 48.06, 4.26, 4.26, 9.57, 9.57, 9.87, 9.87, 4.12, 4.12, 10.08, 10.08, 23.61, 23.61, 9.02, 9.02, 25.89, 25.89, 22.36, 22.36, 40.48, 40.48, 13.79, 13.79, 18.93, 18.93, 30, 30, 20.01, 20.01, 25.19, 25.19, 4.47, 4.47, 9.86, 9.86, 12.93, 12.93, 18.85, 18.85, 48.63, 48.63, 14.37, 14.37, 0.69, 0.69, 24.18, 24.18, 13.8, 13.8, 20.37, 20.37, 22.03, 22.03, 45.46, 45.46, 6.01, 6.01, 9.92, 9.92, 13.35, 13.35, 46.75, 46.75, 26.26, 26.26, 36.98, 36.98, 27.78, 27.78, 11.29, 11.29, 25.13, 25.13, 13.99, 13.99, 48.17, 48.17, 26.12, 26.12, 36.05, 36.05, 14.12, 14.12, 46.47, 46.47, 42.61, 42.61, 1.03, 1.03, 16.1, 16.1, 40.42, 40.42, 28.65, 28.65, 12.42, 12.42, 44.86, 44.86, 0.73, 0.73, 21.39, 21.39, 4.44, 4.44, 39.98, 39.98, 39.01, 39.01, 36.31, 36.31, 5.86, 5.86, 28.89, 28.89, 35.94, 35.94, 7.13, 7.13, 38.66, 38.66, 17.63, 17.63, 44.57, 44.57, 22.77, 22.77, 21.07, 21.07, 33.53, 33.53, 38.45, 38.45, 30.6, 30.6, 8.81, 8.81, 42.4, 42.4, 21.3, 21.3, 1.94, 1.94, 8.65, 8.65, 13.42, 13.42, 45.67, 45.67, 12.59, 12.59, 27.66, 27.66, 30.52, 30.52, 4.23, 4.23, 11.89, 11.89, 6.45, 6.45, 7.68, 7.68, 27.72, 27.72, 40.23, 40.23, 0.06, 0.06, 10.68, 10.68, 26.16, 26.16, 22.87, 22.87, 20, 20, 26.12, 26.12, 38.93, 38.93, 46.5, 46.5, 21.53, 21.53, 28.41, 28.41, 32.87, 32.87, 0.47, 0.47, 8.85, 8.85, 16.59, 16.59, 44.07, 44.07, 29, 29, 19.75, 19.75, 30.68, 30.68, 15.99, 15.99, 4.27, 4.27, 2.22, 2.22, 45.38, 45.38, 3.33, 3.33, 19.34, 19.34, 3.25, 3.25, 9.18, 9.18, 23.95, 23.95, 6.9, 6.9, 16.45, 16.45, 39.17, 39.17, 2.44, 2.44, 26.81, 26.81, 42.52, 42.52, 43.13, 43.13, 1.53, 1.53, 2.17, 2.17, 45.38, 45.38, 9.05, 9.05, 20.36, 20.36, 17.94, 17.94, 18.25, 18.25), + size = c(29, 29, 28, 28, 26, 26, 35, 35, 22, 22, 38, 38, 14, 14, 25, 25, 31, 31, 24, 24, 22, 22, 18, 18, 26, 26, 27, 27, 31, 31, 21, 21, 22, 22, 31, 31, 32, 32, 31, 31, 25, 25, 29, 29, 24, 24, 22, 22, 30, 30, 27, 27, 30, 30, 16, 16, 32, 32, 27, 27, 30, 30, 24, 24, 29, 29, 35, 35, 22, 22, 25, 25, 19, 19, 26, 26, 29, 29, 26, 26, 25, 25, 26, 26, 17, 17, 30, 30, 22, 22, 26, 26, 21, 21, 23, 23, 36, 36, 26, 26, 28, 28, 27, 27, 19, 19, 35, 35, 22, 22, 26, 26, 23, 23, 26, 26, 17, 17, 23, 23, 24, 24, 20, 20, 24, 24, 19, 19, 24, 24, 33, 33, 29, 29, 33, 33, 20, 20, 28, 28, 24, 24, 28, 28, 28, 28, 25, 25, 20, 20, 30, 30, 27, 27, 21, 21, 26, 26, 24, 24, 26, 26, 31, 31, 17, 17, 30, 30, 21, 21, 25, 25, 34, 34, 16, 16, 23, 23, 24, 24, 19, 19, 27, 27, 31, 31, 23, 23, 26, 26, 30, 30, 31, 31, 11, 11, 28, 28, 24, 24, 27, 27, 19, 19, 26, 26, 22, 22, 32, 32, 37, 37, 21, 21, 16, 16, 23, 23, 20, 20, 23, 23, 20, 20, 27, 27, 25, 25, 27, 27, 27, 27, 20, 20, 23, 23, 29, 29, 21, 21, 22, 22, 26, 26, 16, 16, 22, 22, 21, 21, 24, 24, 26, 26, 18, 18, 21, 21, 22, 22, 30, 30, 23, 23, 33, 33, 23, 23, 25, 25, 29, 29, 22, 22, 22, 22, 30, 30, 22, 22, 24, 24, 20, 20, 18, 18, 23, 23, 24, 24, 29, 29, 25, 25, 18, 18, 35, 35, 22, 22, 27, 27, 20, 20, 31, 31, 29, 29, 29, 29, 22, 22, 20, 20, 30, 30, 33, 33, 27, 27, 22, 22, 25, 25, 17, 17, 18, 18, 29, 29, 25, 25, 24, 24, 29, 29, 25, 25, 27, 27, 28, 28, 30, 30, 34, 34, 18, 18, 24, 24, 19, 19, 30, 30, 19, 19, 34, 34, 18, 18, 30, 30, 27, 27, 14, 14, 30, 30, 12, 12, 35, 35), + height = c(166, 166, 159, 159, 172, 172, 171, 171, 163, 163, 173, 173, 162, 162, 163, 163, 174, 174, 161, 161, 158, 158, 157, 157, 159, 159, 160, 160, 162, 162, 156, 156, 159, 159, 164, 164, 166, 166, 172, 172, 168, 168, 171, 171, 162, 162, 168, 168, 159, 159, 165, 165, 159, 159, 166, 166, 156, 156, 164, 164, 157, 157, 168, 168, 161, 161, 167, 167, 166, 166, 162, 162, 160, 160, 159, 159, 159, 159, 167, 167, 161, 161, 161, 161, 157, 157, 169, 169, 174, 174, 168, 168, 167, 167, 168, 168, 165, 165, 160, 160, 169, 169, 167, 167, 163, 163, 172, 172, 164, 164, 168, 168, 155, 155, 171, 171, 173, 173, 159, 159, 172, 172, 164, 164, 170, 170, 173, 173, 163, 163, 163, 163, 162, 162, 153, 153, 166, 166, 167, 167, 167, 167, 159, 159, 158, 158, 168, 168, 171, 171, 163, 163, 170, 170, 159, 159, 165, 165, 162, 162, 168, 168, 165, 165, 164, 164, 165, 165, 163, 163, 171, 171, 157, 157, 154, 154, 163, 163, 158, 158, 166, 166, 166, 166, 163, 163, 167, 167, 169, 169, 168, 168, 157, 157, 171, 171, 165, 165, 166, 166, 162, 162, 170, 170, 169, 169, 161, 161, 160, 160, 179, 179, 163, 163, 162, 162, 156, 156, 168, 168, 165, 165, 163, 163, 165, 165, 159, 159, 160, 160, 159, 159, 155, 155, 157, 157, 163, 163, 159, 159, 167, 167, 169, 169, 166, 166, 162, 162, 160, 160, 159, 159, 168, 168, 167, 167, 161, 161, 172, 172, 174, 174, 171, 171, 156, 156, 158, 158, 168, 168, 166, 166, 173, 173, 161, 161, 156, 156, 167, 167, 161, 161, 164, 164, 162, 162, 160, 160, 163, 163, 163, 163, 167, 167, 160, 160, 162, 162, 169, 169, 170, 170, 152, 152, 168, 168, 164, 164, 166, 166, 164, 164, 161, 161, 165, 165, 162, 162, 160, 160, 166, 166, 164, 164, 172, 172, 171, 171, 171, 171, 166, 166, 170, 170, 159, 159, 161, 161, 164, 164, 176, 176, 164, 164, 166, 166, 170, 170, 170, 170, 168, 168, 158, 158, 160, 160, 161, 161, 168, 168, 165, 165, 163, 163, 168, 168, 173, 173, 165, 165, 161, 161), + sex = c('male', 'male', 'male', 'male', 'male', 'male', 'male', 'male', 'female', 'female', 'female', 'female', 'female', 'female', 'female', 'female', 'female', 'female', 'male', 'male', 'male', 'male', 'male', 'male', 'male', 'male', 'female', 'female', 'female', 'female', 'male', 'male', 'male', 'male', 'female', 'female', 'female', 'female', 'female', 'female', 'male', 'male', 'male', 'male', 'female', 'female', 'male', 'male', 'male', 'male', 'male', 'male', 'male', 'male', 'male', 'male', 'male', 'male', 'male', 'male', 'male', 'male', 'male', 'male', 'male', 'male', 'female', 'female', 'male', 'male', 'female', 'female', 'male', 'male', 'male', 'male', 'female', 'female', 'male', 'male', 'male', 'male', 'female', 'female', 'female', 'female', 'female', 'female', 'male', 'male', 'female', 'female', 'female', 'female', 'male', 'male', 'female', 'female', 'male', 'male', 'male', 'male', 'female', 'female', 'female', 'female', 'female', 'female', 'male', 'male', 'male', 'male', 'male', 'male', 'male', 'male', 'female', 'female', 'male', 'male', 'male', 'male', 'male', 'male', 'female', 'female', 'female', 'female', 'male', 'male', 'male', 'male', 'female', 'female', 'male', 'male', 'male', 'male', 'male', 'male', 'male', 'male', 'female', 'female', 'male', 'male', 'male', 'male', 'male', 'male', 'female', 'female', 'female', 'female', 'male', 'male', 'female', 'female', 'female', 'female', 'male', 'male', 'male', 'male', 'male', 'male', 'male', 'male', 'male', 'male', 'female', 'female', 'female', 'female', 'male', 'male', 'male', 'male', 'female', 'female', 'female', 'female', 'male', 'male', 'male', 'male', 'female', 'female', 'male', 'male', 'female', 'female', 'female', 'female', 'male', 'male', 'female', 'female', 'male', 'male', 'female', 'female', 'male', 'male', 'female', 'female', 'male', 'male', 'male', 'male', 'male', 'male', 'male', 'male', 'male', 'male', 'male', 'male', 'female', 'female', 'female', 'female', 'female', 'female', 'male', 'male', 'male', 'male', 'female', 'female', 'male', 'male', 'female', 'female', 'female', 'female', 'male', 'male', 'male', 'male', 'female', 'female', 'male', 'male', 'female', 'female', 'female', 'female', 'female', 'female', 'male', 'male', 'female', 'female', 'female', 'female', 'male', 'male', 'female', 'female', 'female', 'female', 'male', 'male', 'female', 'female', 'male', 'male', 'male', 'male', 'male', 'male', 'female', 'female', 'female', 'female', 'female', 'female', 'female', 'female', 'male', 'male', 'female', 'female', 'female', 'female', 'male', 'male', 'female', 'female', 'male', 'male', 'male', 'male', 'female', 'female', 'female', 'female', 'male', 'male', 'male', 'male', 'male', 'male', 'female', 'female', 'male', 'male', 'female', 'female', 'female', 'female', 'male', 'male', 'female', 'female', 'female', 'female', 'female', 'female', 'male', 'male', 'female', 'female', 'male', 'male', 'female', 'female', 'female', 'female', 'female', 'female', 'male', 'male', 'female', 'female', 'female', 'female', 'male', 'male', 'female', 'female', 'male', 'male', 'male', 'male', 'male', 'male', 'male', 'male', 'female', 'female', 'male', 'male', 'male', 'male', 'female', 'female', 'male', 'male', 'male', 'male', 'male', 'male', 'female', 'female', 'male', 'male', 'male', 'male', 'male', 'male'), + col = c('red', 'red', 'blue', 'blue', 'red', 'red', 'red', 'red', 'red', 'red', 'blue', 'blue', 'red', 'red', 'red', 'red', 'blue', 'blue', 'red', 'red', 'blue', 'blue', 'red', 'red', 'red', 'red', 'blue', 'blue', 'red', 'red', 'red', 'red', 'red', 'red', 'blue', 'blue', 'red', 'red', 'red', 'red', 'blue', 'blue', 'red', 'red', 'red', 'red', 'blue', 'blue', 'blue', 'blue', 'blue', 'blue', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'blue', 'blue', 'red', 'red', 'blue', 'blue', 'red', 'red', 'red', 'red', 'red', 'red', 'blue', 'blue', 'blue', 'blue', 'red', 'red', 'blue', 'blue', 'red', 'red', 'blue', 'blue', 'red', 'red', 'blue', 'blue', 'red', 'red', 'blue', 'blue', 'blue', 'blue', 'blue', 'blue', 'blue', 'blue', 'blue', 'blue', 'red', 'red', 'blue', 'blue', 'red', 'red', 'blue', 'blue', 'blue', 'blue', 'red', 'red', 'blue', 'blue', 'red', 'red', 'red', 'red', 'blue', 'blue', 'blue', 'blue', 'blue', 'blue', 'red', 'red', 'red', 'red', 'blue', 'blue', 'red', 'red', 'blue', 'blue', 'blue', 'blue', 'red', 'red', 'blue', 'blue', 'blue', 'blue', 'red', 'red', 'red', 'red', 'blue', 'blue', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'blue', 'blue', 'blue', 'blue', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'blue', 'blue', 'red', 'red', 'blue', 'blue', 'blue', 'blue', 'red', 'red', 'blue', 'blue', 'red', 'red', 'red', 'red', 'red', 'red', 'blue', 'blue', 'blue', 'blue', 'blue', 'blue', 'blue', 'blue', 'blue', 'blue', 'red', 'red', 'red', 'red', 'blue', 'blue', 'blue', 'blue', 'blue', 'blue', 'red', 'red', 'blue', 'blue', 'blue', 'blue', 'red', 'red', 'blue', 'blue', 'blue', 'blue', 'red', 'red', 'red', 'red', 'red', 'red', 'blue', 'blue', 'blue', 'blue', 'red', 'red', 'blue', 'blue', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'blue', 'blue', 'blue', 'blue', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'blue', 'blue', 'blue', 'blue', 'blue', 'blue', 'red', 'red', 'blue', 'blue', 'blue', 'blue', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'blue', 'blue', 'blue', 'blue', 'blue', 'blue', 'red', 'red', 'red', 'red', 'red', 'red', 'blue', 'blue', 'blue', 'blue', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'blue', 'blue', 'red', 'red', 'red', 'red', 'blue', 'blue', 'blue', 'blue', 'red', 'red', 'blue', 'blue', 'blue', 'blue', 'blue', 'blue', 'blue', 'blue', 'blue', 'blue', 'blue', 'blue', 'blue', 'blue', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'red', 'blue', 'blue', 'blue', 'blue', 'blue', 'blue', 'red', 'red', 'red', 'red')) + +# Add distend distbegin to the data +tot.data <- Distance:::create_bins(tot.data, c(0,10,25,50)) + +# Code to simulate data above +# library(dsims) +# +# covs <- list() +# covs$size <- list(list(distribution = "poisson", lambda = 25)) +# covs$height <- list(list(distribution = "normal", mean = 165, sd = 5)) +# covs$sex <- data.frame(level = c("male", "female"), +# prob = c(0.5, 0.5)) +# covs$col <- data.frame(level = c("blue", "red"), +# prob = c(0.4, 0.6)) +# +# # Define the population description (this time using the density to determine +# # the population size) +# popdesc <- make.population.description(covariates = covs, +# N = 250) +# +# cov.param <- list() +# cov.param$size <- c(log(1.02)) +# cov.param$height <- c(log(1.005)) +# cov.param$sex <- data.frame(level = c("male", "female"), +# param = c(log(1.5), 0)) +# cov.param$col <- data.frame(level = c("blue", "red"), +# param = c(log(1.5), 0)) +# +# +# # define the detecability +# detect <- make.detectability(key.function = "hn", +# scale.param = 5, +# cov.param = cov.param, +# truncation = 50) +# +# plot(detect, popdesc) +# +# simulation <- make.simulation(reps = 1, +# population.description = popdesc, +# detectability = detect) +# +# # run an example survey to check the setup +# survey <- run.survey(simulation) +# +# survey2 <- run.survey(survey, make.region()) +# +# ddata1 <- survey@dist.data +# ddata2 <- survey2@dist.data +# +# head(ddata1) +# head(ddata2) +# +# # Add observer and detected columns +# ddata1$observer <- 1 +# ddata2$observer <- 2 +# ddata1$detected <- 1 +# ddata2$detected <- 1 +# +# # Combine the two datasets +# tot.data <- rbind(ddata1, ddata2) +# +# # Check which objects only have 1 entry (all need 2) +# index <- which(table(tot.data$individual) == 1) +# index <- as.numeric(names(index)) +# +# # Add rows for animals that were detected by one observer but not the other +# tmp <- tot.data[tot.data$individual %in% index,] +# +# # These animals were only detected by one observer +# tmp$detected <- 0 +# new.observer <- tmp$observer +# new.observer <- ifelse(new.observer==1,2,1) +# tmp$observer <- new.observer +# +# tot.data <- rbind(tot.data, tmp) +# +# # Order the data - by observer +# index <- order(tot.data$observer) +# tot.data <- tot.data[index,] +# +# # Order by individual +# index <- order(tot.data$individual) +# tot.data <- tot.data[index,] +# +# # Remove object column +# tot.data <- tot.data[,-1] +# # Rename individual column to object +# names(tot.data)[1] <- "object" + +test_that("gof when no df",{ + + fit.mr <- ddf(data = tot.data, + method = "io", + dsmodel = ~mcds(key="hr", formula = ~size+height+sex), + mrmodel = ~glm(formula = ~size+height+sex+col), + meta.data = list(binned = TRUE, breaks = c(0,10,25,50))) + + # Check it knows there are insufficient df for pooled chisquare + tmp <- ddf.gof(fit.mr) + expect_true(is.na(tmp$chisquare$pooled.chi$df)) + expect_true(is.na(tmp$chisquare$pooled.chi$p)) + +}) From e8ac8d0802fe41b4830b1f22b26f095c00c222e7 Mon Sep 17 00:00:00 2001 From: Laura Marshall Date: Mon, 3 Feb 2025 11:32:12 +0000 Subject: [PATCH 02/14] fix other MR gof total chi sq df and p values reference #96 --- R/gof.io.fi.R | 17 +++++++++++++---- R/gof.rem.R | 15 ++++++++++++--- R/gof.rem.fi.R | 15 ++++++++++++--- R/gof.trial.R | 15 ++++++++++++--- R/gof.trial.fi.R | 15 ++++++++++++--- R/print.ddf.gof.R | 1 + 6 files changed, 62 insertions(+), 16 deletions(-) diff --git a/R/gof.io.fi.R b/R/gof.io.fi.R index 857aa7e5..e24db7c1 100644 --- a/R/gof.io.fi.R +++ b/R/gof.io.fi.R @@ -73,7 +73,17 @@ gof.io.fi <- function(model,breaks=NULL,nc=NULL){ p.1 <- NA df.2 <- NA p.2 <- NA - + + # Calculate the pooled chi-square + df.pool <- 3*nc-length(model$par)-1 + + if(df.pool <= 0){ + df.pool <- NA + p.pool <- NA + }else{ + p.pool <- 1-pchisq(chisq.1+chisq.2, df.pool) + } + return(list(chi1=list(observed=observed.count.1, expected=expected.1, chisq=chisq.1, @@ -85,7 +95,6 @@ gof.io.fi <- function(model,breaks=NULL,nc=NULL){ p=p.2, df=df.2), pooled.chi=list(chisq=chisq.1+chisq.2, - df=3*nc-length(model$par)-1, - p=1-pchisq(chisq.1+chisq.2, - 3*nc-length(model$par)-1)))) + df=df.pool, + p=p.pool))) } diff --git a/R/gof.rem.R b/R/gof.rem.R index 993caec7..4c5b7190 100644 --- a/R/gof.rem.R +++ b/R/gof.rem.R @@ -77,6 +77,16 @@ gof.rem <- function(model,breaks=NULL,nc=NULL){ }else{ p.2 <- 1-pchisq(chisq.2,df.2) } + + # Calculate the pooled chi-square + df.pool <- 2*nc-length(model$par)-1 + + if(df.pool <= 0){ + df.pool <- NA + p.pool <- NA + }else{ + p.pool <- 1-pchisq(chisq.1+chisq.2, df.pool) + } return(list(chi1=list(observed=observed.count.1, expected=expected.1, @@ -89,7 +99,6 @@ gof.rem <- function(model,breaks=NULL,nc=NULL){ p=p.2, df=df.2), pooled.chi=list(chisq=chisq.1+chisq.2, - df=2*nc-length(model$par)-1, - p=1-pchisq(chisq.1+chisq.2, - 2*nc-length(model$par)-1)))) + df=df.pool, + p=p.pool))) } diff --git a/R/gof.rem.fi.R b/R/gof.rem.fi.R index f39c66c6..785474c1 100644 --- a/R/gof.rem.fi.R +++ b/R/gof.rem.fi.R @@ -69,6 +69,16 @@ gof.rem.fi <- function(model,breaks=NULL,nc=NULL){ p.1 <- NA df.2 <- NA p.2 <- NA + + # Calculate the pooled chi-square + df.pool <- 2*nc-length(model$par)-1 + + if(df.pool <= 0){ + df.pool <- NA + p.pool <- NA + }else{ + p.pool <- 1-pchisq(chisq.1+chisq.2, df.pool) + } return(list(chi1=list(observed=observed.count.1, expected=expected.1, @@ -81,7 +91,6 @@ gof.rem.fi <- function(model,breaks=NULL,nc=NULL){ p=p.2, df=df.2), pooled.chi=list(chisq=chisq.1+chisq.2, - df=2*nc-length(model$par)-1, - p=1-pchisq(chisq.1+chisq.2, - 2*nc-length(model$par)-1)))) + df=df.pool, + p=p.pool))) } diff --git a/R/gof.trial.R b/R/gof.trial.R index 83914199..9a0e92c7 100644 --- a/R/gof.trial.R +++ b/R/gof.trial.R @@ -81,6 +81,16 @@ gof.trial <- function(model, breaks=NULL, nc=NULL){ }else{ p.2 <- 1-pchisq(chisq.2,df.2) } + + # Calculate the pooled chi-square + df.pool <- 2*nc-length(model$par)-1 + + if(df.pool <= 0){ + df.pool <- NA + p.pool <- NA + }else{ + p.pool <- 1-pchisq(chisq.1+chisq.2, df.pool) + } return(list(chi1=list(observed=observed.count.1, expected=expected.1, @@ -93,7 +103,6 @@ gof.trial <- function(model, breaks=NULL, nc=NULL){ p=p.2, df=df.2), pooled.chi=list(chisq=chisq.1+chisq.2, - df=2*nc-length(model$par)-1, - p=1-pchisq(chisq.1+chisq.2, - 2*nc-length(model$par)-1)))) + df=df.pool, + p=p.pool))) } diff --git a/R/gof.trial.fi.R b/R/gof.trial.fi.R index 7f715657..9a2806c7 100644 --- a/R/gof.trial.fi.R +++ b/R/gof.trial.fi.R @@ -69,6 +69,16 @@ gof.trial.fi <- function(model,breaks=NULL,nc=NULL){ p.1 <- NA df.2 <- NA p.2 <- NA + + # Calculate the pooled chi-square + df.pool <- 2*nc-length(model$par)-1 + + if(df.pool <= 0){ + df.pool <- NA + p.pool <- NA + }else{ + p.pool <- 1-pchisq(chisq.1+chisq.2, df.pool) + } return(list(chi1=list(observed=observed.count.1, expected=expected.1, @@ -81,7 +91,6 @@ gof.trial.fi <- function(model,breaks=NULL,nc=NULL){ p=p.2, df=df.2), pooled.chi=list(chisq=chisq.1+chisq.2, - df=2*nc-length(model$par)-1, - p=1-pchisq(chisq.1+chisq.2, - 2*nc-length(model$par)-1)))) + df=df.pool, + p=p.pool))) } diff --git a/R/print.ddf.gof.R b/R/print.ddf.gof.R index 4bbaaf7d..f17f0dd2 100644 --- a/R/print.ddf.gof.R +++ b/R/print.ddf.gof.R @@ -30,6 +30,7 @@ print.ddf.gof <- function(x, digits=3, ...){ if(!is.null(gof$chisquare)){ cat("\nChi-square tests\n") + # This is NULL when it is a single observer model if(!is.null(gof$chisquare$chi2)){ cat("\nDistance sampling component:\n") } From b83c2b42bc8f3b84a1ea2f71a521d6d7036a8b94 Mon Sep 17 00:00:00 2001 From: Laura Marshall Date: Tue, 18 Mar 2025 11:49:46 +0000 Subject: [PATCH 03/14] Update NEWS.md --- NEWS.md | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/NEWS.md b/NEWS.md index d18495a8..3703c3ad 100644 --- a/NEWS.md +++ b/NEWS.md @@ -3,7 +3,7 @@ Bug Fixes * Fixed formatting issue in flnl.grad help -* Fixed output for gof test when degrees of freedom for MR io models is <= 0 +* Fixed output for gof test when degrees of freedom for MR models is <= 0 # mrds 3.0.0 From b548f4e34f29b9da877e3ae908146053846693cb Mon Sep 17 00:00:00 2001 From: Laura Marshall Date: Thu, 19 Jun 2025 12:21:42 +0100 Subject: [PATCH 04/14] Documentation update Reference issue [#186](https://github.com/DistanceDevelopment/Distance/issues/186#issuecomment-2496236254) --- R/check.mono.R | 2 +- man/check.mono.Rd | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/R/check.mono.R b/R/check.mono.R index 054afb4a..c3297efa 100644 --- a/R/check.mono.R +++ b/R/check.mono.R @@ -21,7 +21,7 @@ #' @param strict if \code{TRUE} (default) the detection function must be #' "strictly" monotone, that is that (\code{g(x[i])<=g(x[i-1])}) over the whole #' range (left to right truncation points). -#' @param n.pts number of equally-spaced points between left and right +#' @param n.pts number of points between left and right #' truncation at which to evaluate the detection function (default 100) #' @param tolerance numerical tolerance for monotonicity checks (default 1e-8) #' @param plot plot a diagnostic highlighting the non-monotonic areas (default diff --git a/man/check.mono.Rd b/man/check.mono.Rd index ef1d83be..11d0e73f 100644 --- a/man/check.mono.Rd +++ b/man/check.mono.Rd @@ -20,7 +20,7 @@ check.mono( "strictly" monotone, that is that (\code{g(x[i])<=g(x[i-1])}) over the whole range (left to right truncation points).} -\item{n.pts}{number of equally-spaced points between left and right +\item{n.pts}{number of points between left and right truncation at which to evaluate the detection function (default 100)} \item{tolerance}{numerical tolerance for monotonicity checks (default 1e-8)} From 77ed1aba2b7c4f5d1d31c394891d45f46b4222b5 Mon Sep 17 00:00:00 2001 From: Laura Marshall Date: Tue, 1 Jul 2025 20:36:43 +0100 Subject: [PATCH 05/14] No longer plot distance points for binned analyses Reference [Distance Issue 144](https://github.com/DistanceDevelopment/Distance/issues/144) --- DESCRIPTION | 2 +- NEWS.md | 1 + R/plot.ds.R | 3 ++- 3 files changed, 4 insertions(+), 2 deletions(-) diff --git a/DESCRIPTION b/DESCRIPTION index 306d29e3..a2f3a6be 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -18,7 +18,7 @@ Description: Animal abundance estimation via conventional, multiple covariate fitting is performed via maximum likelihood. Also included are diagnostics and plotting for fitted detection functions. Abundance estimation is via a Horvitz-Thompson-like estimator. -Version: 3.0.0.9006 +Version: 3.0.0.9007 URL: https://github.com/DistanceDevelopment/mrds/ BugReports: https://github.com/DistanceDevelopment/mrds/issues Depends: diff --git a/NEWS.md b/NEWS.md index 5ad05e58..ec024386 100644 --- a/NEWS.md +++ b/NEWS.md @@ -7,6 +7,7 @@ Bug Fixes * Now displays a warning if the user tries to fit a detection function with covariates using MCDS.exe which is not either a half-normal or a hazard rate model. (Issue #113) * Fixed so that the MCDS.exe does not try to fit a negative exponential in place of a gamme key function. (Issue #113) * Now issues warnings when there is only a single transect and varflag option is 1 or 2. (Issue #115) +* No longer plots data points on detection functions for binned analyses, only shows the interval bins. # mrds 3.0.0 diff --git a/R/plot.ds.R b/R/plot.ds.R index 9a0d9197..80223d21 100644 --- a/R/plot.ds.R +++ b/R/plot.ds.R @@ -382,7 +382,8 @@ plot.ds <- function(x, which=2, breaks=NULL, nc=NULL, ldots$y <- linevalues do.call(lines, ldots) - if(showpoints){ + # Don't show point if the data are binned + if(showpoints && !binned){ jitter.p <- rnorm(length(point_vals), 1, jitter.v[1]) pdots <- dots pdots$x <- xmat$distance From 8589899f502180e5687aaa821cb2de15c18b8978 Mon Sep 17 00:00:00 2001 From: Laura Marshall Date: Thu, 3 Jul 2025 13:41:16 +0100 Subject: [PATCH 06/14] Update dht.se documentation Closes #117 --- R/dht.se.R | 17 +++++++++++++++-- man/dht.se.Rd | 17 +++++++++++++++-- 2 files changed, 30 insertions(+), 4 deletions(-) diff --git a/R/dht.se.R b/R/dht.se.R index 819830ee..6812596c 100644 --- a/R/dht.se.R +++ b/R/dht.se.R @@ -9,7 +9,7 @@ #' \itemize{ #' \item variation due to uncertainty from estimation of the detection #' function parameters; -#' \item variation in abundance due to random sample selection; +#' \item variation in abundance due to random sample selection. #' } #' The first component (model parameter uncertainty) is computed using a delta #' method estimate of variance (\insertCite{huggins1989;nobrackets}{mrds}; \insertCite{huggins1991;nobrackets}{mrds}; \insertCite{borchers1998;nobrackets}{mrds}) in @@ -34,7 +34,7 @@ #' is used in the detection function; #' \item \code{2} is the default and uses the encounter rate estimator #' \eqn{\hat{N}/L} (estimated abundance per unit transect) suggested by -#' \insertCite{innes2002;textual}{mrds} and \insertCite{marques2004;textual}{mrds} +#' \insertCite{innes2002;textual}{mrds} and \insertCite{marques2004;textual}{mrds}. #' } #' #' In general if any covariates are used in the models, the default @@ -66,6 +66,19 @@ #' of the assumed log-normal distribution because the components are #' multiplicative. For combining df for the sampling variance across regions #' they are weighted by the variance because it is a sum across regions. +#' +#' The coefficient of variation (CV) associated with the abundance estimates is calculated based on the following formula for the \code{varflag} options 1 and 2: +#' +#' \code{varflag=1} +#' +#' \deqn{CV(\hat{N}) = \sqrt{\left(\frac{\sqrt{n}}{n}\right)^2+CV(\hat{p})^2}} +#' +#' \code{varflag=2} +#' +#' \deqn{CV(\hat{N}) = \sqrt{\left(\frac{\sqrt{\hat{N}}}{\hat{N}}\right)^2+CV(\hat{p})^2}} +#' where n is the number of observations, \eqn{\hat{N}} is the estimated +#' abundance and \eqn{\hat{p}} is the average probability of detection for +#' an animal in the covered area. #' #' A non-zero correlation between regional estimates can occur from using a #' common detection function across regions. This is reflected in the diff --git a/man/dht.se.Rd b/man/dht.se.Rd index c0e885d3..8d59fbf1 100644 --- a/man/dht.se.Rd +++ b/man/dht.se.Rd @@ -47,7 +47,7 @@ The variance has two components: \itemize{ \item variation due to uncertainty from estimation of the detection function parameters; - \item variation in abundance due to random sample selection; + \item variation in abundance due to random sample selection. } The first component (model parameter uncertainty) is computed using a delta method estimate of variance (\insertCite{huggins1989;nobrackets}{mrds}; \insertCite{huggins1991;nobrackets}{mrds}; \insertCite{borchers1998;nobrackets}{mrds}) in @@ -72,7 +72,7 @@ to calculate encounter rate: is used in the detection function; \item \code{2} is the default and uses the encounter rate estimator \eqn{\hat{N}/L} (estimated abundance per unit transect) suggested by - \insertCite{innes2002;textual}{mrds} and \insertCite{marques2004;textual}{mrds} + \insertCite{innes2002;textual}{mrds} and \insertCite{marques2004;textual}{mrds}. } In general if any covariates are used in the models, the default @@ -105,6 +105,19 @@ of the assumed log-normal distribution because the components are multiplicative. For combining df for the sampling variance across regions they are weighted by the variance because it is a sum across regions. +The coefficient of variation (CV) associated with the abundance estimates is calculated based on the following formula for the \code{varflag} options 1 and 2: + +\code{varflag=1} + +\deqn{CV(\hat{N}) = \sqrt{\left(\frac{\sqrt{n}}{n}\right)^2+CV(\hat{p})^2}} + +\code{varflag=2} + +\deqn{CV(\hat{N}) = \sqrt{\left(\frac{\sqrt{\hat{N}}}{\hat{N}}\right)^2+CV(\hat{p})^2}} +where n is the number of observations, \eqn{\hat{N}} is the estimated +abundance and \eqn{\hat{p}} is the average probability of detection for +an animal in the covered area. + A non-zero correlation between regional estimates can occur from using a common detection function across regions. This is reflected in the correlation matrix of the regional and total estimates which is given in the From 508cd9b42707a75c8f70bea91be64a2a7e99f627 Mon Sep 17 00:00:00 2001 From: Laura Marshall Date: Thu, 3 Jul 2025 13:44:18 +0100 Subject: [PATCH 07/14] Fix a plotting issue I introduced an issue when turning off plotting the points for binned data. Reference mrds issue [144](https://github.com/DistanceDevelopment/Distance/issues/144) --- R/plot.ds.R | 12 +++++++----- 1 file changed, 7 insertions(+), 5 deletions(-) diff --git a/R/plot.ds.R b/R/plot.ds.R index 80223d21..ecc5d8e1 100644 --- a/R/plot.ds.R +++ b/R/plot.ds.R @@ -175,13 +175,15 @@ plot.ds <- function(x, which=2, breaks=NULL, nc=NULL, zgrid <- matrix(rep(z[1,], length(xgrid)), byrow=TRUE, ncol=sum(zdim)) } + # Check if it is a binned analysis + if(is.null(model$meta.data$binned)){ + binned <- FALSE + }else{ + binned <- model$meta.data$binned + } + # create intervals of distance (breaks) for the chosen number of classes (nc). if(is.null(breaks)){ - if(is.null(model$meta.data$binned)){ - binned <- FALSE - }else{ - binned <- model$meta.data$binned - } if(binned){ breaks <- model$ds$aux$breaks nc <- length(breaks)-1 From f23de2981cb41684474e0c61c0a15aeddba3131f Mon Sep 17 00:00:00 2001 From: Laura Marshall Date: Thu, 3 Jul 2025 13:45:22 +0100 Subject: [PATCH 08/14] Fix rem model gof degrees of freedom Reference #96 --- R/gof.rem.R | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/R/gof.rem.R b/R/gof.rem.R index 4c5b7190..d7bcc888 100644 --- a/R/gof.rem.R +++ b/R/gof.rem.R @@ -70,7 +70,7 @@ gof.rem <- function(model,breaks=NULL,nc=NULL){ p.1 <- 1-pchisq(chisq.1,df.1) } - df.2 <- 2*nc-length(model$mr$par) + df.2 <- nc-length(model$mr$par) if(df.2<=0){ df.2 <- NA p.2 <- NA From 1cf470b4fa9faf71192f111051918d751b8f103b Mon Sep 17 00:00:00 2001 From: Laura Marshall Date: Thu, 3 Jul 2025 14:08:17 +0100 Subject: [PATCH 09/14] Modify the bysample table Reference issue #60 and #100 --- R/dht.R | 23 +++++++++++++++++------ man/dht.Rd | 8 ++++++-- 2 files changed, 23 insertions(+), 8 deletions(-) diff --git a/R/dht.R b/R/dht.R index ada18b89..593633af 100644 --- a/R/dht.R +++ b/R/dht.R @@ -152,8 +152,12 @@ #' #' The list structure of clusters and individuals are the same: #' \item{bysample}{\code{data.frame} giving results for each sample; -#' \code{Nchat} is the estimated abundance within the sample and \code{Nhat} is -#' scaled by surveyed area/covered area within that region} +#' \code{Sample.Area} is the covered area associated with the sampler, +#' \code{n} is the number of detections on the sampler, +#' \code{Nhat} is the estimated abundance within the sample, and +#' \code{Dhat} is \eqn{\frac{Nhat}{\sum{Sample.Area}}} so that summing +#' these values gives the overall density estimates.} +#' #' \item{summary}{\code{data.frame} of summary statistics for each region and #' total} #' \item{N}{\code{data.frame} of estimates of abundance for each region and @@ -223,15 +227,22 @@ dht <- function(model, region.table, sample.table, obs.table=NULL, subset=NULL, bysample.table <- with(Nhat.by.sample, data.frame(Region = Region.Label, + Area = Area, Sample = Sample.Label, Effort = Effort.x, Sample.Area = CoveredArea, - Area = Area, n = n, - Nhat = Nhat, - Nchat = Nhat*CoveredArea/Area)) + Nhat = Nhat*CoveredArea/Area)) - bysample.table$Dhat <- bysample.table$Nchat/bysample.table$Sample.Area + # Calculate density contributions (these can be summed to give overall density) + bysample.table$Dhat <- bysample.table$Nhat/bysample.table$Sample.Area + # Now update CoveredArea so it's per sample + if(point){ + bysample.table$Sample.Area <- pi*bysample.table$Effort*width^2 - pi*bysample.table$Effort*left^2 + }else{ + bysample.table$Sample.Area <- 2*bysample.table$Effort*(width-left) + } + Nhat.by.region <- as.numeric(by(Nhat.by.sample$Nhat, Nhat.by.sample$Region.Label, sum)) diff --git a/man/dht.Rd b/man/dht.Rd index 66bf9cd7..7187ef0c 100644 --- a/man/dht.Rd +++ b/man/dht.Rd @@ -50,8 +50,12 @@ list object of class \code{dht} with elements: The list structure of clusters and individuals are the same: \item{bysample}{\code{data.frame} giving results for each sample; -\code{Nchat} is the estimated abundance within the sample and \code{Nhat} is -scaled by surveyed area/covered area within that region} +\code{Sample.Area} is the covered area associated with the sampler, +\code{n} is the number of detections on the sampler, +\code{Nhat} is the estimated abundance within the sample, and +\code{Dhat} is \eqn{\frac{Nhat}{\sum{Sample.Area}}} so that summing +these values gives the overall density estimates.} + \item{summary}{\code{data.frame} of summary statistics for each region and total} \item{N}{\code{data.frame} of estimates of abundance for each region and From b26a245b5d041a7873c6c2ae6c5c5131dafc02d1 Mon Sep 17 00:00:00 2001 From: Laura Marshall Date: Sat, 5 Jul 2025 02:20:59 +0100 Subject: [PATCH 10/14] Additional documentation for dht summary tables Reference [Distance Issue #182](https://github.com/DistanceDevelopment/Distance/issues/182) --- R/dht.R | 4 +++- man/dht.Rd | 4 +++- 2 files changed, 6 insertions(+), 2 deletions(-) diff --git a/R/dht.R b/R/dht.R index 593633af..901f7a5f 100644 --- a/R/dht.R +++ b/R/dht.R @@ -159,7 +159,9 @@ #' these values gives the overall density estimates.} #' #' \item{summary}{\code{data.frame} of summary statistics for each region and -#' total} +#' total. Note that the summary statistics give a general summary of the data +#' and may use more basic calculations than those used in the abundance +#' and density calculations.} #' \item{N}{\code{data.frame} of estimates of abundance for each region and #' total} #' \item{D}{\code{data.frame} of estimates of density for each region and total} diff --git a/man/dht.Rd b/man/dht.Rd index 7187ef0c..bc667cc8 100644 --- a/man/dht.Rd +++ b/man/dht.Rd @@ -57,7 +57,9 @@ The list structure of clusters and individuals are the same: these values gives the overall density estimates.} \item{summary}{\code{data.frame} of summary statistics for each region and -total} +total. Note that the summary statistics give a general summary of the data +and may use more basic calculations than those used in the abundance +and density calculations.} \item{N}{\code{data.frame} of estimates of abundance for each region and total} \item{D}{\code{data.frame} of estimates of density for each region and total} From 9ea5479641000e20ca9ec7d321115a80e22699e8 Mon Sep 17 00:00:00 2001 From: Laura Marshall Date: Sat, 5 Jul 2025 02:22:30 +0100 Subject: [PATCH 11/14] Add check on object field of obs.table Reference [Distance Issue #165](https://github.com/DistanceDevelopment/Distance/issues/165) --- R/dht.R | 9 +++++++++ 1 file changed, 9 insertions(+) diff --git a/R/dht.R b/R/dht.R index 901f7a5f..ded462e4 100644 --- a/R/dht.R +++ b/R/dht.R @@ -419,6 +419,15 @@ dht <- function(model, region.table, sample.table, obs.table=NULL, subset=NULL, ervar=ifelse(model$meta.data$point, "P2", "R2"), areas.supplied=FALSE) + + # Input checking + # Check that the object field is numeric in the obs.table + if(!is.null(obs.table)){ + if(!is.numeric(obs.table$object)){ + warning("Please ensure the object field in the obs.table is numeric, cannot perform dht calculations.", call. = FALSE) + return(NULL) + } + } # jll 18-Nov-04; the following allows for a subset statement to be added to # create obs.table from model data rather than creating obs.table separately. From 40021cb637203427abcbd6e9c713c11caf5a4821 Mon Sep 17 00:00:00 2001 From: Laura Marshall Date: Sat, 5 Jul 2025 02:22:43 +0100 Subject: [PATCH 12/14] Update NEWS and DESCRIPTION --- DESCRIPTION | 2 +- NEWS.md | 10 +++++++++- 2 files changed, 10 insertions(+), 2 deletions(-) diff --git a/DESCRIPTION b/DESCRIPTION index a2f3a6be..6da16797 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -18,7 +18,7 @@ Description: Animal abundance estimation via conventional, multiple covariate fitting is performed via maximum likelihood. Also included are diagnostics and plotting for fitted detection functions. Abundance estimation is via a Horvitz-Thompson-like estimator. -Version: 3.0.0.9007 +Version: 3.0.1 URL: https://github.com/DistanceDevelopment/mrds/ BugReports: https://github.com/DistanceDevelopment/mrds/issues Depends: diff --git a/NEWS.md b/NEWS.md index ec024386..ccbe3a1c 100644 --- a/NEWS.md +++ b/NEWS.md @@ -3,11 +3,19 @@ Bug Fixes * Fixed formatting issue in flnl.grad help -* Fixed output for gof test when degrees of freedom for MR models is <= 0 +* Fixed output for gof test when degrees of freedom for MR models is <= 0 (Issue #96) * Now displays a warning if the user tries to fit a detection function with covariates using MCDS.exe which is not either a half-normal or a hazard rate model. (Issue #113) * Fixed so that the MCDS.exe does not try to fit a negative exponential in place of a gamme key function. (Issue #113) * Now issues warnings when there is only a single transect and varflag option is 1 or 2. (Issue #115) * No longer plots data points on detection functions for binned analyses, only shows the interval bins. +* Add a check and a warning when the object field of the obs.table is not numeric. (Distance Issue #165) + + +Enhancements + +* Documentation of dht.se to show variance calculations for varflag options (Issue #117) +* Modified the bysample table returned from dht (Issues #60 and #100) +* Clarified documentation regarding summary table output from dht. # mrds 3.0.0 From 1d2e5af9d2fcb0c91541e89065be5a000bedaeb5 Mon Sep 17 00:00:00 2001 From: Laura Marshall Date: Sat, 5 Jul 2025 02:33:51 +0100 Subject: [PATCH 13/14] Fix spelling errors --- NEWS.md | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/NEWS.md b/NEWS.md index ccbe3a1c..8e5487f4 100644 --- a/NEWS.md +++ b/NEWS.md @@ -166,7 +166,7 @@ Bug Fixes * rescaling parameters were not correct, now fixed. Thanks to Laura Marshall for spotting this. * coefficients are called coefficients (not a mixture of coefficients and parameters) in summary() results -* speed-up in io.fi models (thanks to Winston Chang's profvis, showing many unecessary calls to model.matrix) +* speed-up in io.fi models (thanks to Winston Chang's profvis, showing many unnecessary calls to model.matrix) * plot.ds now has a pdf= option to plot the probability density function (for point transect models only) * assign.par, create.ddfobj and detfct are now exported, so it can be used by dsm (though shouldn't be used by anything else!) * fixed bug in left truncation where probability of detection was not calculated correctly. Thanks to Jason Roberts for pointing this out! @@ -227,7 +227,7 @@ BUG FIXES BUG FIXES -* Standardisation was being applied to detection functions (such that g(0)=1) when there were no adjustments (which is uneccesary) but also caused issues when using gamma detection functions as this should be calculated at g(apex) instead. Standardisation code has been removed for when there are no adjustments and the correct scaling used for the gamma when there are. Thanks to Thomas Doniol-Valcroze for alerting us to this bug. +* Standardisation was being applied to detection functions (such that g(0)=1) when there were no adjustments (which is unnecessary) but also caused issues when using gamma detection functions as this should be calculated at g(apex) instead. Standardisation code has been removed for when there are no adjustments and the correct scaling used for the gamma when there are. Thanks to Thomas Doniol-Valcroze for alerting us to this bug. * Partial name-matching in dht was fixed. Produced warning but not error. NEW FEATURES @@ -264,7 +264,7 @@ BUG FIXES * uniform+cosine detection functions were ignored when using monotonicity constraints, now they can be used together * mono.strict=TRUE didn't automatically turn on mono=TRUE, extra logic to correct this -* montonicity constraints did not use standardised (g(x)/g(0) detection functions, so if g(x)>1 monotonicity constraints were voilated. Now standardised detection functions are used. Thanks to Len Thomas for noticing this bug. +* montonicity constraints did not use standardised (g(x)/g(0) detection functions, so if g(x)>1 monotonicity constraints were violated. Now standardised detection functions are used. Thanks to Len Thomas for noticing this bug. # mrds 2.1.4-3 @@ -275,7 +275,7 @@ BUG FIX CHANGES * general documentation updates -* simplication and re-structuring of internals +* simplification and re-structuring of internals # mrds 2.1.4-3 From a52acb7c949b7a40a3faf1a1fa03b27b95820ee5 Mon Sep 17 00:00:00 2001 From: Laura Marshall Date: Sat, 5 Jul 2025 03:10:55 +0100 Subject: [PATCH 14/14] pkgdown build --- docs/404.html | 2 +- docs/LICENSE-text.html | 2 +- docs/articles/index.html | 2 +- docs/articles/mrds-golftees.html | 4 ++-- docs/authors.html | 6 +++--- docs/index.html | 2 +- docs/news/index.html | 17 ++++++++++++----- docs/pkgdown.yml | 6 +++--- docs/reference/AIC.ddf.html | 2 +- docs/reference/DeltaMethod.html | 2 +- docs/reference/NCovered.html | 2 +- docs/reference/add.df.covar.line.html | 2 +- docs/reference/adj.check.order.html | 2 +- docs/reference/adj.cos.html | 2 +- docs/reference/adj.herm.html | 2 +- docs/reference/adj.poly.html | 2 +- docs/reference/adj.series.grad.cos.html | 2 +- docs/reference/adj.series.grad.herm.html | 2 +- docs/reference/adj.series.grad.poly.html | 2 +- docs/reference/apex.gamma.html | 2 +- docs/reference/assign.default.values.html | 2 +- docs/reference/assign.par.html | 2 +- docs/reference/average.line.cond.html | 2 +- docs/reference/average.line.html | 2 +- docs/reference/book.tee.data.html | 2 +- docs/reference/calc.se.Np.html | 2 +- docs/reference/cdf.ds.html | 2 +- docs/reference/cds.html | 2 +- docs/reference/check.bounds.html | 2 +- docs/reference/check.mono.html | 4 ++-- docs/reference/coef.ds.html | 2 +- docs/reference/compute.Nht.html | 2 +- docs/reference/covered.region.dht.html | 2 +- docs/reference/create.bins.html | 2 +- docs/reference/create.command.file.html | 2 +- docs/reference/create.ddfobj.html | 2 +- docs/reference/create.model.frame.html | 2 +- docs/reference/create.varstructure.html | 2 +- docs/reference/ddf.ds.html | 2 +- docs/reference/ddf.gof.html | 2 +- docs/reference/ddf.html | 2 +- docs/reference/ddf.io.fi.html | 2 +- docs/reference/ddf.io.html | 2 +- docs/reference/ddf.rem.fi.html | 2 +- docs/reference/ddf.rem.html | 2 +- docs/reference/ddf.trial.fi.html | 2 +- docs/reference/ddf.trial.html | 2 +- docs/reference/det.tables.html | 2 +- docs/reference/detfct.fit.html | 2 +- docs/reference/detfct.fit.opt.html | 2 +- docs/reference/dht.deriv.html | 2 +- docs/reference/dht.html | 14 ++++++++++---- docs/reference/dht.se.html | 15 ++++++++++++--- docs/reference/distpdf.grad.html | 2 +- docs/reference/distpdf.html | 2 +- docs/reference/ds.function.html | 2 +- docs/reference/flnl.constr.grad.neg.html | 2 +- docs/reference/flnl.grad.html | 2 +- docs/reference/flnl.html | 2 +- docs/reference/flt.var.html | 2 +- docs/reference/g0.html | 2 +- docs/reference/getpar.html | 2 +- docs/reference/gof.ds.html | 2 +- docs/reference/gstdint.html | 2 +- docs/reference/histline.html | 2 +- docs/reference/index.html | 2 +- docs/reference/integratedetfct.logistic.html | 2 +- docs/reference/integratelogistic.analytic.html | 2 +- docs/reference/integratepdf.grad.html | 2 +- docs/reference/integratepdf.html | 2 +- docs/reference/io.glm.html | 2 +- docs/reference/is.linear.logistic.html | 2 +- docs/reference/is.logistic.constant.html | 2 +- docs/reference/keyfct.grad.hn.html | 2 +- docs/reference/keyfct.grad.hz.html | 2 +- docs/reference/keyfct.th1.html | 2 +- docs/reference/keyfct.th2.html | 2 +- docs/reference/keyfct.tpn.html | 2 +- docs/reference/lfbcvi.html | 2 +- docs/reference/lfgcwa.html | 2 +- docs/reference/logLik.ddf.html | 2 +- docs/reference/logisticbyx.html | 2 +- docs/reference/logisticbyz.html | 2 +- docs/reference/logisticdetfct.html | 2 +- docs/reference/logisticdupbyx.html | 2 +- docs/reference/logisticdupbyx_fast.html | 2 +- docs/reference/logit.html | 2 +- docs/reference/mcds.html | 2 +- docs/reference/mcds_dot_exe.html | 2 +- docs/reference/mrds-package.html | 2 +- docs/reference/mrds_opt.html | 2 +- docs/reference/nlminb_wrapper.html | 2 +- docs/reference/p.det.html | 2 +- docs/reference/p.dist.table.html | 2 +- docs/reference/parse.optimx.html | 2 +- docs/reference/pdot.dsr.integrate.logistic.html | 2 +- docs/reference/plot.det.tables.html | 2 +- docs/reference/plot.ds.html | 2 +- docs/reference/plot.io.fi.html | 2 +- docs/reference/plot.io.html | 2 +- docs/reference/plot.rem.fi.html | 2 +- docs/reference/plot.rem.html | 2 +- docs/reference/plot.trial.fi.html | 2 +- docs/reference/plot.trial.html | 2 +- docs/reference/plot_cond.html | 2 +- docs/reference/plot_layout.html | 2 +- docs/reference/plot_uncond.html | 2 +- docs/reference/predict.ds.html | 2 +- docs/reference/print.ddf.gof.html | 2 +- docs/reference/print.ddf.html | 2 +- docs/reference/print.det.tables.html | 2 +- docs/reference/print.dht.html | 2 +- docs/reference/print.p_dist_table.html | 2 +- docs/reference/print.summary.ds.html | 2 +- docs/reference/print.summary.io.fi.html | 2 +- docs/reference/print.summary.io.html | 2 +- docs/reference/print.summary.rem.fi.html | 2 +- docs/reference/print.summary.rem.html | 2 +- docs/reference/print.summary.trial.fi.html | 2 +- docs/reference/print.summary.trial.html | 2 +- docs/reference/prob.deriv.html | 2 +- docs/reference/prob.se.html | 2 +- docs/reference/process.data.html | 2 +- docs/reference/pronghorn.html | 2 +- docs/reference/ptdata.distance.html | 2 +- docs/reference/ptdata.dual.html | 2 +- docs/reference/ptdata.removal.html | 2 +- docs/reference/ptdata.single.html | 2 +- docs/reference/qqplot.ddf.html | 2 +- docs/reference/rem.glm.html | 2 +- docs/reference/rescale_pars.html | 2 +- docs/reference/sample_ddf.html | 2 +- docs/reference/setbounds.html | 2 +- docs/reference/setcov.html | 2 +- docs/reference/setinitial.ds.html | 2 +- docs/reference/sim.mix.html | 2 +- docs/reference/solvecov.html | 2 +- docs/reference/stake77.html | 2 +- docs/reference/stake78.html | 2 +- docs/reference/summary.ds.html | 2 +- docs/reference/summary.io.fi.html | 2 +- docs/reference/summary.io.html | 2 +- docs/reference/summary.rem.fi.html | 2 +- docs/reference/summary.rem.html | 2 +- docs/reference/summary.trial.fi.html | 2 +- docs/reference/summary.trial.html | 2 +- docs/reference/survey.region.dht.html | 2 +- docs/reference/test.breaks.html | 2 +- docs/reference/varn.html | 2 +- docs/search.json | 2 +- 150 files changed, 187 insertions(+), 165 deletions(-) diff --git a/docs/404.html b/docs/404.html index ea4b0fee..7264b32e 100644 --- a/docs/404.html +++ b/docs/404.html @@ -28,7 +28,7 @@ mrds - 3.0.0.9006 + 3.0.1