Why some iSSA's estimated coefficients are presented as NA?

[AlineGiroux]

I am using amt to run iSSA and I am not sure why some coefficients of the model are presented as NA in the model's summary.
I am using the very same model with different species (different data) and this happen to some of them (some species have no NA and present estimated values for all model coefficients).
Can you help me please? @bsmity13 @jmsigner

[bsmity13]

It looks like your variable map is a factor with levels forest and open, right? The NA estimates all appear in interactions with either c1 or tpt, which I'm guessing are continuous. And it also seems like those NA only occur for the forest category.

R typically uses "dummy" coding to represent factor variables, with k-1 coefficients to represent a factor with k levels, with the 1st level captured by the intercept. The thing with these conditional logistic regression models is that there is no intercept, so I think R (or really, the survival package) tends to create too many levels for factor interactions. I don't fully understand why that is.

Essentially, it seems like forest should be the reference category, and so the model is overparameterized with the interaction terms. In other words, the variable map is a logical variable for map == "open". So when map == "forest", the result is 0, and 0 times whatever the values of c1 or tpt are is still 0, so you have no variation to work with.

You can fix all of this by making the interaction variables yourself. So instead of map, create a variable that is map_open with map == "open", then use that in interaction with tpt and c1. You might want to consider main effects of tpt and c1 in your model formula, as well. The main effects will represent the selection for tpt and c1 when map == "forest", and the interaction slopes will add together with the main effects to represent selection when map == "open".

Have I missed anything? Does that make sense? Does it fix the problem?

[AlineGiroux]

Thanks a lot Brian! It makes a lot of sense.
Yes, map is a factor with levels forest and open. I tried use open as reference running ssf$map <- relevel(factor(ssf$map),"open") before run the model. But it seems like the model is using forest as reference anyway. How can I change that?

Yes, the variables c1 and tpt are time of day (continuous) and environmental temperature (continuous). I didn't consider main effects of c1 and tpt because this variables do not vary within a step (they vary in time, but not in space). This is why they only appear in interaction with other variables.

This is the model:

model <- fit_issf(case_~
                 #BIASED KERNEL HABITAT
                 map +
                 map:c1 +
                 map:tpt +
                 map:c1:tpt +
                 
                 #MOVEMENT FREE KERNEL
                 sl_ + log_sl_ +
                 sl_:c1 + log_sl_:c1 +
                 sl_:tpt + log_sl_:tpt +
                 sl_:map + log_sl_:map +
                 sl_:c1:tpt + log_sl_:c1:tpt +
                 sl_:c1:map + log_sl_:c1:map +
                 sl_:tpt:map + log_sl_:tpt:map +
                 sl_:c1:tpt:map + log_sl_:c1:tpt:map +
                 
                 cos_ta_ +
                 cos_ta_:c1 + 
                 cos_ta_:tpt +
                 cos_ta_:map +
                 cos_ta_:c1:tpt + 
                 cos_ta_:c1:map +
                 cos_ta_:tpt:map +
                 cos_ta_:c1:tpt:map +
                 
                 strata(step_id_),data=ssf,model=T)

I am not sure if I understood how to include the 'map' in the model. Using open as reference, should I create a column in ssf with 0 for open and 1 for forest? And then include this column in the model instead of map? And if I had more than two categories in map?

[AlineGiroux]

I tried to follow your suggestion.

ssf$map_forest <- ssf$map == "forest"

model1 <- fit_issf(case_~
                 #BIASED KERNEL HABITAT
                 map +
                 c1 + map_forest:c1 +
                 tpt + map_forest:tpt +
                 c1:tpt + map_forest:c1:tpt +
                 
                 #MOVEMENT FREE KERNEL
                 sl_ + log_sl_ +
                 sl_:c1 + log_sl_:c1 +
                 sl_:tpt + log_sl_:tpt +
                 sl_:map_forest + log_sl_:map_forest +
                 sl_:c1:tpt + log_sl_:c1:tpt +
                 sl_:c1:map_forest + log_sl_:c1:map_forest +
                 sl_:tpt:map_forest + log_sl_:tpt:map_forest +
                 sl_:c1:tpt:map_forest + log_sl_:c1:tpt:map_forest +
                 
                 cos_ta_ +
                 cos_ta_:c1 + 
                 cos_ta_:tpt +
                 cos_ta_:map_forest +
                 cos_ta_:c1:tpt + 
                 cos_ta_:c1:map_forest +
                 cos_ta_:tpt:map_forest +
                 cos_ta_:c1:tpt:map_forest +
                 
                 strata(step_id_),data=ssf,model=T)

But then I got NA to the main effects of c1 and tpt (as expected because they do not vary within steps). Is that right? And if I had more than two categories?

[jmsigner]

A couple of things you could try:

1. `ssf$map_forest <- as.nummeric(ssf$map == "forest")``
2. Change

#BIASED KERNEL HABITAT
                 map +

to

#BIASED KERNEL HABITAT
                 map_forest +

3. Use cland tpt only in interactions, as they do not vary per stratum.

[AlineGiroux]

If I remove the main effects of c1 and tpt from the model, then I get NA for the same interactions that I was getting in the very beginning.

model2 <- fit_issf(case_~
                 #BIASED KERNEL HABITAT
                 map +
                 map_forest:c1 +
                 map_forest:tpt +
                 map_forest:c1:tpt +
                 
                 #MOVEMENT FREE KERNEL
                 sl_ + log_sl_ +
                 sl_:c1 + log_sl_:c1 +
                 sl_:tpt + log_sl_:tpt +
                 sl_:map_forest + log_sl_:map_forest +
                 sl_:c1:tpt + log_sl_:c1:tpt +
                 sl_:c1:map_forest + log_sl_:c1:map_forest +
                 sl_:tpt:map_forest + log_sl_:tpt:map_forest +
                 sl_:c1:tpt:map_forest + log_sl_:c1:tpt:map_forest +
                 
                 cos_ta_ +
                 cos_ta_:c1 + 
                 cos_ta_:tpt +
                 cos_ta_:map_forest +
                 cos_ta_:c1:tpt + 
                 cos_ta_:c1:map_forest +
                 cos_ta_:tpt:map_forest +
                 cos_ta_:c1:tpt:map_forest +
                 
                 strata(step_id_),data=ssf,model=T)

Then, comparing this model with the previous one, we can notice that the effect of map_forestTRUE:c1 is equal to map_forestFALSE:c1 *(-1). The same logic can be observed to map_forestTRUE:tpt and map_forestTRUE:c1:tpt. It makes sense for two categories, but if I had more then two? Should I keep c1 and tpt in the model as main effects even knowing that they do not vary within steps?

Thank you very very much!

[jmsigner]

Did you also try this ssf$map_forest <- as.nummeric(ssf$map == "forest")?

[AlineGiroux]

Sorry, I didn't. I only saw your message after sending mine. Will do now.

[jmsigner]

Should I keep c1 and tpt in the model as main effects even knowing that they do not vary within steps?

No

[AlineGiroux]

Perfect!
using ssf$map_forest <- as.numeric(ssf$map == "forest"), and running:

model3 <- fit_issf(case_~
                 #BIASED KERNEL HABITAT
                 map_forest +
                 map_forest:c1 +
                 map_forest:tpt +
                 map_forest:c1:tpt +
                 
                 #MOVEMENT FREE KERNEL
                 sl_ + log_sl_ +
                 sl_:c1 + log_sl_:c1 +
                 sl_:tpt + log_sl_:tpt +
                 sl_:map_forest + log_sl_:map_forest +
                 sl_:c1:tpt + log_sl_:c1:tpt +
                 sl_:c1:map_forest + log_sl_:c1:map_forest +
                 sl_:tpt:map_forest + log_sl_:tpt:map_forest +
                 sl_:c1:tpt:map_forest + log_sl_:c1:tpt:map_forest +
                 
                 cos_ta_ +
                 cos_ta_:c1 + 
                 cos_ta_:tpt +
                 cos_ta_:map_forest +
                 cos_ta_:c1:tpt + 
                 cos_ta_:c1:map_forest +
                 cos_ta_:tpt:map_forest +
                 cos_ta_:c1:tpt:map_forest +
                 
                 strata(step_id_),data=ssf,model=T)

summary(model3)

I got

Thanks @bsmity13 and @jmsigner!

Any idea how I could do if I had more than two categories on the map?

[jmsigner]

Any idea how I could do if I had more than two categories on the map?

You can do the dummy encoding for several categories by hand and/or use model.matrix(): A small example here:

dat <- data.frame(
  landcover = c("forest", "forest", "water", "water", "other")
)

# manually
landcover.other <- as.numeric(dat$landcover == "other")
landcover.water <- as.numeric(dat$landcover == "water")



# Use model matrix
dat1 <- model.matrix(~ landcover, data = dat) |> as.data.frame()
dat1[, -1]

[AlineGiroux]

Understood. Many thanks!