library("blm") data("aarp") fit <- blm(bladder70~female * smoke_status, aarp, weigh=aarp$w) coef(fit)*1000 summary(fit) AIC(fit) which.at.boundary(fit) confint(fit)*1000 model.formula(fit) all.vars(model.formula(fit)) risk.types <- unique(subset(aarp, select=all.vars(model.formula(fit)))) risk.types <- subset(risk.types, bladder70==0) risk.types predict(fit, risk.types, se=TRUE)*1000 Rsquared(fit) risk <- crude.risk(bladder70~fiber.centered, data=aarp, weights=aarp$w) head(risk) risk.exposure.plot(object=risk, scale = 1000, las=1, col = "royalblue", pch = 19, ylab = "Crude risk (per 1,000)", xlab = "Avg. Fiber Consumption (Centered)") # USE LEXPIT TO EXPAND MODEL formula.linear <- bladder70~female * smoke_status formula.expit <- bladder70~redmeat+fiber.centered+I(fiber.centered^2) fit <- lexpit(formula.linear, formula.expit, aarp, weight=aarp$w) summary(fit) AIC(fit) # TEST FORM fit.both <- lexpit(update(formula.linear,~.+fiber.centered+I(fiber.centered^2)), formula.expit, data=aarp, weight=aarp$w) summary(fit.both) # IMPROVED FIT Rsquared(fit) EO(fit, aarp$educ) gof(fit) CIs <- confint(fit) CIs[1:7,]*1000 CIs[8:11,] expit(CIs[8,])*10000