problem with estimating N using region.N

[Dina Matyukhina]

Dear group members,

I'm relatively new to SECR. And I get overwhelmed every time an error message pops up and I can't find a solution. So I'd really appreciate your suggestions and thoughts. 

I'm analyzing six-year Amur leopard data. I fit a number of models with spatial covariates on density and with additive or interaction terms of the session (year). 

The top model did not include a year covariate. 

model9b_ppo120<-secr.fit(input, model=list(D~sd1000, g0~h2,sigma~h2,pmix~1), method="Nelder-Mead", detectfn=0, mask=mask2, hcov = "Sex", sessioncov=session_cov, start=model9a_ppo120)

When I tried to apply region.N command 

region.N(model9b_ppo120, region=mask2, session=1)

 It gave me the error

Error in pmixn[x, ] <- ifelse(knownclass > 1, ifelse(knownclass == (x +  :

  number of items to replace is not a multiple of replacement length

And it gives me the same error for each session except session 5

region.N(model9b_ppo120, region=mask2, session=5)

         estimate    SE.estimate       lcl                     ucl               n

E.N  83.83885    4.795464          74.95447        93.7763      77

R.N 104.65917    NaN                  NaN                 NaN            77

I assume this model should produce the mean N estimate, and therefore for each year, E.N would be the same. The only thing that puzzles me is what special about session 5 might be or what might be wrong with the rest of the data given this type of error.

To estimate N for each year separately I fit the next model

model26_ppo120<-secr.fit(input, model=list(D~year, g0~h2,sigma~h2,pmix~1), method="Nelder-Mead", detectfn=0, mask=mask2, hcov = "Sex", sessioncov=session_cov)

But it gave me the following error for every single session, including session 5

region.N(model26b_ppo120, region=mask2, session=5)

Error in `contrasts<-`(`*tmp*`, value = contr.funs[1 + isOF[nn]]) :

contrasts can be applied only to factors with 2 or more levels

I use R version i 4.0.1. and secr version 4.2.2.

Thank you!

Regards,

Dina

[Murray Efford]

Dina

I understand this all can be aggravating! At first I thought you might just be pushing a limited dataset too far, but I was wrong. It seems there is a bug in region.N for multi-session models with non-null hcov (I can reproduce the problem with the ovenCHp dataset). I'll have to investigate further - give me a couple of days. (Also, I don't think specifying pmix~1 has any effect - just leave that out of your model specification).

Murray

[Murray Efford]

Sorry this took longer than I expected. I tracked down the bug and fixed it in the pre-release version of secr 4.3.2 that is available from

https://www.otago.ac.nz/density/zip/secr_4.3.2.zip (windows binary) or https://www.otago.ac.nz/density/zip/secr_4.3.2.tar.gz (source).

Murray

[Dina Matyukhina]

Hello Murray,

Thank you a lot for taking the time to fix it and sending the link.

Unfortunately, I got exactly the same error after installing the new version of the package.

library(secr)
This is secr 4.3.2 pre-release. For overview type ?secr

I specified sessioncov for my six-year data as

session_cov=data.frame(constant=c("1", "1", "1", "1", "1","1"),
                       year=c("1", "2", "3", "4", "5","6"),
                       trend=c(1,2,3,4,5,6))

The model is 

model25b<-secr.fit(input, model=list(D~year, g0~1,sigma~1), method="Nelder-Mead", detectfn=0,
                   mask=mask2, sessioncov=session_cov, start=model25a)

But when I tried to apply region.N to this model, I got the error

region.N(model25b, region=mask2, session=3)

Error in `contrasts<-`(`*tmp*`, value = contr.funs[1 + isOF[nn]]) :
  contrasts can be applied only to factors with 2 or more levels

So any model with a year covariate produces this type of error.

I'm just wondering, whether there is something wrong with my dataset. However, when I first ran the analysis a couple of years ago using the same code it was fine. 

Regards,

Dina

среда, 4 ноября 2020 г. в 10:39:45 UTC+10, murray...@gmail.com:

[Murray Efford]

I'm sorry I overlooked the second, unrelated problem in your original post.

I think this can be blamed on the changed default in R >=4.0 for the data.frame argument 'stringsAsFactors'. I predict it will work if you use

session_cov = data.frame(constant = c("1", "1", "1", "1", "1","1"),
                       year = factor(c("1", "2", "3", "4", "5","6")),
                       trend = c(1,2,3,4,5,6))

You can even substitute this retrospectively in a fitted model:

model25$sessioncov <- session_cov

region.N( model2)

Why you are specifying an alternative 'year' covariate when there is a 1:1 relation between years and sessions? D~session has the same effect.

Murray

[Dina Matyukhina]

Thank you, Murray, I really appreciate it. That worked! And the first problem I mentioned in the initial post was solved after fixing the bug. Sorry about the confusion. 

I somehow thought that I need to specify a year covariate in this way. Now I know I shouldn't. 

But in the case I have two seasons within one session (year) can I specify sessions like this?

session_cov=data.frame(year=c("1", "1", "2", "2","3", "3", "4", "4"),

                       season=c("1", "2", "1", "2","1", "2", "1", "2"))

Regards,

Dina

пятница, 6 ноября 2020 г. в 05:48:47 UTC+10, murray...@gmail.com:

[Murray Efford]

Sure - that's what session covariates are good for.

For the record, I've fixed the code for density models so character-valued session covariates should automatically be converted to factors. Now in the pre-release version on the website as given before.

Murray