The hardware and bandwidth for this mirror is donated by dogado GmbH, the Webhosting and Full Service-Cloud Provider. Check out our Wordpress Tutorial.
If you wish to report a bug, or if you are interested in having us mirror your free-software or open-source project, please feel free to contact us at mirror[@]dogado.de.

BayesGP: Fitting Model with Partial Likelihood

library(BayesGP)

Partial Likelihood Models in BayesGP

There are currently two implemented models in BayesGP that use the partial likelihood function for inference: the case-crossover model and the Cox Proportional Hazard (Coxph) model.

Case-Crossover Model

With BayesGP, one can specify the argument family to "cc", "casecrossover" or "CaseCrossover" to fit a case-crossover model.

Here we will use a simulated dataset:

data <- as.data.frame(ccData)
data$exposure <- data$exposure
mod <- model_fit(formula = case ~ f(x = exposure, 
                                    model = "IWP", 
                                    order = 2, k = 30,
                                    initial_location = median(data$exposure), 
                                    sd.prior = list(prior = "exp", param = list(u = 1, alpha = 0.5), h = 1)),
                 family = "cc",
                 strata = "subject",
                 weight = NULL,
                 data = data,
                 method = "aghq")

To take a look at its result:

true_effect <- function(x) {3 *(x^2 - .5^2)}
plot(mod)

lines(I(true_effect(seq(0,1,by = 0.1)) - true_effect(median(data$exposure))) ~ seq(0,1,by = 0.1), col = "red")

Here the true effect used to simulate the data is shown as the red line. It is important to know that for case-crossover model, the intercept parameter and the strata level effects will not be identifiable.

CoxPH Model

For Cox Proportional Hazard Model, one can specify the argument family to "coxph to fit a CoxPH model with its partial likelihood.

Here we will illustrate with the kidney example from the survival package.

data <- survival::kidney
head(data)
#>   id time status age sex disease frail
#> 1  1    8      1  28   1   Other   2.3
#> 2  1   16      1  28   1   Other   2.3
#> 3  2   23      1  48   2      GN   1.9
#> 4  2   13      0  48   2      GN   1.9
#> 5  3   22      1  32   1   Other   1.2
#> 6  3   28      1  32   1   Other   1.2
mod <- model_fit(formula = time ~ age + sex + f(x = id, 
                                    model = "IID", 
                                    sd.prior = list(prior = "exp", param = list(u = 1, alpha = 0.5))),
                 family = "coxph",
                 cens = "status",
                 data = data,
                 method = "aghq")

Take a look at the posterior for each fixed effect:

samps_age <- sample_fixed_effect(mod, variables = "age")
samps_sex <- sample_fixed_effect(mod, variables = "sex")
par(mfrow = c(1,2))
hist(samps_age, main = "Samples for effect of age", xlab = "Effect")
hist(samps_sex, main = "Samples for effect of sex", xlab = "Effect")

par(oldpar)

These binaries (installable software) and packages are in development.
They may not be fully stable and should be used with caution. We make no claims about them.
Health stats visible at Monitor.