NAs produced by resDirM

[mtvincent-NOAA]

I've been using your package to calculate the one step ahead residuals for age and length compositions that are fit by a Dirichlet-multinomial distribution. I noticed that the redDirM function does not work if you have any zeros in the observations. Therefore, I had added a small constant to the observed and predicted to try to overcome this issue. However, I still occasionally get the error message: In qnorm(pred$Fx - U * pred$px) : NaNs produced
Do you know what causes this or how it can be avoided?
I can send you some examples if necessary.

[vtrijoulet]

Hi, a short reproducible example would be useful.
When using the Dirichlet multinomial, observations can be zeros but not the alphas. Also the observations need to be integers.

[mtvincent-NOAA]

okay. I didn't realize that the alphas needed to be integers. What are your suggestions if the predictions (alphas?) have zeros in the matrix?

[mtvincent-NOAA]

Also I seem to get errors with a zero in the observation with this example.

a<-matrix(rep(1000*c(0,.2,.1,.1,.1,.3), 100), nrow=6)
o<-rdirM(100,1000, a[,1])
res<-resDirM(o,a)

Warning messages:
1: In nlminb(obj$par, obj$fn, obj$gr) : NA/NaN function evaluation
2: In nlminb(obj$par, obj$fn, obj$gr) : NA/NaN function evaluation

[calbertsen]

The issue in your example is not the zeros in the observations (o), but the row of zeros in the concentration parameters / alphas (a).
For example,

a<-matrix(rep(1000*c(1e-5,.2,.1,.1,.1,.3), 100), nrow=6)
o<-rdirM(100,1000, a[,1])
res<-resDirM(o,a)

Still has zeros in the observations, but does not give a warning.

[mtvincent-NOAA]

I finally had enough time to put together a reproducible example. See the example below where the rows are the ages and the columns are the years. I converted the proportions observed to a number by multiplying the sample size and rounding to an integer. The pred is the model predicted proportions at age by year, that mostly sum to 1.
obs= matrix(c( 0, 0, 0, 0, 0, 0, 1,
44, 24, 18, 79, 79, 60, 36,
521, 576, 221, 272, 523, 464, 194,
634, 989, 642, 504, 448, 501, 226,
642, 602, 379, 655, 422, 199, 146,
327, 377, 172, 207, 260, 153, 59,
56, 115, 63, 63, 65, 38, 7,
8, 18, 19, 37, 29, 14, 3,
1, 4, 1, 8, 14, 2, 0,
0, 0, 0, 0, 2, 1, 1,
0, 0, 1, 0, 0, 1, 0,
0, 0, 0, 0, 0, 0, 0),ncol=7,byrow=TRUE)

pred= matrix(c(2.243346e-03, 1.909236e-03, 2.363994e-03, 0.0027393200, 0.0022767180, 2.062993e-03, 0.0082168130,
3.472664e-02, 2.776374e-02, 2.445347e-02, 0.0340818100, 0.0360470700, 2.904686e-02, 0.0838180600,
2.805715e-01, 3.034186e-01, 2.497479e-01, 0.2397383000, 0.3120862000, 3.107834e-01, 0.2558258000,
2.978304e-01, 3.150176e-01, 3.517982e-01, 0.3047442000, 0.2707289000, 3.318746e-01, 0.2855923000,
2.216589e-01, 1.727425e-01, 1.920385e-01, 0.2276582000, 0.1732007000, 1.521490e-01, 0.1943218000,
1.281959e-01, 1.201921e-01, 1.012401e-01, 0.1176402000, 0.1277495000, 9.334646e-02, 0.1090646000,
2.741251e-02, 5.051172e-02, 5.717102e-02, 0.0464461200, 0.0543231400, 5.471914e-02, 0.0390211500,
5.704165e-03, 6.897629e-03, 1.842212e-02, 0.0201531500, 0.0156509900, 1.862296e-02, 0.0152938600,
1.260815e-03, 1.208329e-03, 2.269673e-03, 0.0059317660, 0.0059820950, 4.951991e-03, 0.0062582910,
2.928516e-04, 2.574657e-04, 3.867069e-04, 0.0007120237, 0.0017070510, 1.846384e-03, 0.0018668180,
7.042939e-05, 5.869743e-05, 8.094871e-05, 0.0001192019, 0.0002012374, 5.177739e-04, 0.0004741795,
3.248053e-05, 2.240695e-05, 2.738082e-05, 0.0000356042, 0.0000463281, 7.848201e-05, 0.0002463507),ncol=7,byrow=TRUE)

z1=compResidual::resDirM(obs,pred)

Warning messages:
1: In nlminb(obj$par, obj$fn, obj$gr) : NA/NaN function evaluation
2: In nlminb(obj$par, obj$fn, obj$gr) : NA/NaN function evaluation

To assist with the debugging below is my sessionInfo()
R version 4.1.1 (2021-08-10)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows 10 x64 (build 19045)

Matrix products: default

locale:
[1] LC_COLLATE=English_United States.1252
[2] LC_CTYPE=English_United States.1252
[3] LC_MONETARY=English_United States.1252
[4] LC_NUMERIC=C
[5] LC_TIME=English_United States.1252

attached base packages:
[1] stats graphics grDevices utils datasets methods base

other attached packages:
[1] FishGraph_2.0.2 tseries_0.10-53 nortest_1.0-4 lmtest_0.9-40
[5] zoo_1.8-11 compResidual_0.0.1

loaded via a namespace (and not attached):
[1] compiler_4.1.1 Matrix_1.5-3 tools_4.1.1 TMB_1.9.2
[5] xts_0.13.0 curl_5.0.0 quantmod_0.4.20 grid_4.1.1
[9] TTR_0.24.3 lattice_0.20-44 quadprog_1.5-8

[vtrijoulet]

Hi,

It is true that the resDirM function can be unstable in extreme cases, this was corrected in commit #0a1e28c.
Please reinstall the package following the instructions on the main page.

Double-check you are giving the correct parameters to the function because your alphas sum to one (so likely the ps in Trijoulet et al. 2023 rather than the alphas).

Cheers,

Vanessa