Complex effects (interactions, polynomials, etc.) in writeJagsModel.xxx() about hydenet HOT 11 CLOSED

I think we may have a solution to this in the broom package released to CRAN last month. It will take a bit of rewriting on my part, but Robinson managed to do something I've been struggling to figure out for years.

from hydenet.

I’m back.

Looks appealing. I like the example in their vignette using nonlinear least squares. If we go this route, we may eventually be interested in implementing this as well.

I thought I’d fit a polynomial/interaction model and see what it looked like:

require(HydeNet)
require(broom)

data(PE)

g <- glm(pe ~ poly(d.dimer,2)*angio, family="binomial", data=PE)
tidy(g)

from hydenet.

Now, if I understand this correctly, the formula I would want to implement is

pe ~ -4.12 + 268.92*d.dimer - 0.83 * d.dimer * 3.05 * [angio == 2] + 22.17 * d.dimer * [angio==2] + 2.07 * d.dimer * [angio==2]

Is that correct?

from hydenet.

You didn't have the squares for the two terms involving d.dimer^2

pe ~ -4.12 + 268.92*d.dimer - 0.83 * (d.dimer^2) + 3.05 * [angio == 2] + 22.17 * d.dimer * [angio==2] + 2.07 * (d.dimer^2) * [angio==2]

from hydenet.

The code below works.

library(HydeNet)
g1 <- lm(wells ~ 1, data=PE)
g2 <- glm(pe ~ wells, data=PE, family="binomial")
g3 <- lm(d.dimer ~ pe + pregnant, data=PE)
g4 <- xtabs(~ pregnant, data=PE)
g5 <- glm(angio ~ pe, data=PE, family="binomial")
g6 <- glm(treat ~ d.dimer + angio, data=PE, family="binomial")
g7 <- glm(death ~ pe + treat, data=PE, family="binomial")

bagOfModels <- list(g1,g2,g3,g4,g5,g6,g7)

bagNet <- HydeNetwork(bagOfModels)
bagNet <- setNode(bagNet,
                  treat,
                  nodeFormula = treat ~ poly(d.dimer, 2) + angio,
                  nodeFitter='glm',
                  fitterArgs = list(family=binomial),
                  p = fromData(),
                  fitModel=TRUE)
writeNetworkModel(bagNet, pretty=TRUE)

I've generated two bugs, however. When I use fitModel=FALSE, writeNetworkModel fails.

The second bug is strange. if I use print(bagNet), it prints just fine. but if I just type bagNet, it crashes the system. That concerns me a bit.

Also, we'll have to be careful about passing model objects, it seems. At the moment I don't know that I have a good way to tell that d.dimer is the same node as poly(d.dimer, 2). If you pass the formula in my setNode above in the original model object, it creates two nodes for d.dimer. In theory, it wouldn't be hard to strip out functions, I suppose, but how far down that rabbit hole do we want to go? Do we want to accommodate nested functions? will we take just any function? What if it is a function of two nodes? I'm going to have to think on that one a little bit.

from hydenet.

One more note on this: I haven't incorporated polynomial terms in multinomial logistic regressions. Mostly because there isn't a tidy method for these models. I suspect I'll have to write that method before I can do the polynomial terms in HydeNet.

from hydenet.

Note also that your code above involves a main-effects model, not one that included the interaction term poly(d.dimer,2) : angio

from hydenet.

I suggest we cover the rabbit hole with mulch and move on to decision & utility nodes. That's what I would do at home anyways.

from hydenet.

library(HydeNet)
g1 <- lm(wells ~ 1, data=PE)
g2 <- glm(pe ~ wells, data=PE, family="binomial")
g3 <- lm(d.dimer ~ pe + pregnant, data=PE)
g4 <- xtabs(~ pregnant, data=PE)
g5 <- glm(angio ~ pe, data=PE, family="binomial")
g6 <- glm(treat ~ d.dimer + angio, data=PE, family="binomial")
g7 <- glm(death ~ pe + treat, data=PE, family="binomial")

bagOfModels <- list(g1,g2,g3,g4,g5,g6,g7)

bagNet <- HydeNetwork(bagOfModels)
bagNet <- setNode(bagNet,
                  treat,
                  nodeFormula = treat ~ poly(d.dimer, 2) * angio,
                  p = fromData())

writeNetworkModel(bagNet, pretty=TRUE)

It would have worked with the interaction if I remembered when I was using strings and when I was using regular expressions.

from hydenet.

Reopening this issue, to be addressed after v1.0. And below is another example, with which I don't know what to do.

library(HydeNet)
net <- HydeNetwork(~ card1.ace | card1
                   + card2.ace | card2
                   + initialPoints | card1*card2*card1.ace*card2.ace
                   + hit1 | initialPoints*dealerUpcard
                   + card3 | hit1
                   + card3.ace | card3
                   + pointsAfterCard3 | initialPoints*card3*card3.ace
                   + hit2 | pointsAfterCard3*dealerUpcard
                   + card4 | hit2
                   + card4.ace | card4
                   + pointsAfterCard4 | pointsAfterCard3*card4*card4.ace
                   + hit3 | pointsAfterCard4*dealerUpcard
                   + card5 | hit3
                   + card5.ace | card5
                   + pointsAfterCard5 | pointsAfterCard4*card5*card5.ace
                   + playerFinalPoints | initialPoints*hit1*pointsAfterCard3
                   *hit2*pointsAfterCard4*hit3*pointsAfterCard5
                   + dealerFinalPoints | dealerUpcard
                   + payoff | playerFinalPoints*dealerFinalPoints)

data(BlackJackTrain)

glm.hit1 <- glm(hit1 ~ initialPoints + I(dealerUpcard %in% c("9","10-K","A")),
                data = BlackJackTrain, family="binomial")
glm.hit2 <- glm(hit2 ~ pointsAfterCard3 + I(dealerUpcard %in% c("9","10-K","A")),
                data = BlackJackTrain, family="binomial")
glm.hit3 <- glm(hit3 ~ pointsAfterCard4 + I(dealerUpcard %in% c("9","10-K","A")),
                data = BlackJackTrain, family="binomial")

net <- setNodeModels(net, glm.hit1, glm.hit2, glm.hit3)

Error in parse(text = x) : <text>:1:77: unexpected numeric constant
1: hit1 ~ ilogit(8.16617   + 1.93073 * I(dealerUpcard %in% c("9", "10-K", "A"))TRUE

from hydenet.

This is caused be HydeNet not recognizing the node I(dealerUpcard %in% c("9", "10-k", "A")). This was one of my big fears about translating the formulas is that it can be a bit difficult to strip out the node name from within a function passed to a model. Eventually, I think I can reach a point where if I assume exactly one node name in a function, I can do it effectively.

But that's what is happening here.

from hydenet.