hojsgaard / grim Goto Github PK

Great package Soren! And here comes the first issue:

I am mainly interested in the specification and parameters of the joint probability function of the best (gRim::stepwise.iModel) log-linear model. (I need the join probability function p(i) for later use.)

It seems impossible to get the parameters of the log-linear model from dModel returned by dmod()? As far as I see, dmod() calls fit.dModel() which in turn calls loglin().

loglin() has an argument param = FALSE. Would it be possible to set param = TRUE in fit.dModel() without breaking stuff?

Or am I seeing it wrong and is there a better way to extract the model specification?

Edit:

Parameters can be extracted from best model m.opt with the following code:
loglin(table = m.opt$datainfo$data, margin = m.opt$glist, param = TRUE)$param

Obviously it introduces extra computation, but works.

Now I am looking for a way to define log(p(i)). Would gRbase be the correct place for such a function? Input would be a dModel object.

It seems the current implementation of gRim::stepwise.iModel() only supports starting model-selection from a full main effects model gRim::dmod(~.^1, data = data) if one wants to include all the variables in the selection procedure.

Is there any theoretical reason for this? In other words, couldn't we just as well start from log p(i) = mu and let AIC do the variable selection?

Reason I ask is that starting from a full main effects model limits model selection to low dimensional contingency tables. (e.g. table() will not work on data.frames with more than 31 binary columns since the amount of cells explodes to 2^31.)