Regression Coefficient Extraction from qrjoint Model Fit

Share:

Description

Post process MCMC output from qrjoint to create summaries of intercept and slope function estimates

Usage

1
2
3
4
 
## S3 method for class 'qrjoint'
coef(object, burn.perc = 0.5, nmc = 200, plot = FALSE, show.intercept = TRUE, 
        reduce = TRUE, ...)

Arguments

object

a fitted model of the class qrjoint.

burn.perc

a positive fraction indicating what fraction of the saved draws are to be discarded as burn-in

nmc

integer giving the number of samples, post burn-in, to be used in Monte Carlo averaging

plot

logical indicating if plots are to be made

show.intercept

whether to plot the intercept curve when plot = TRUE

reduce

logical indicating if the tail-expanded grid of tau values is to be reduced to the regular increment grid

...

limited plotting parameters that are passed onto the call of getBands

Value

Extracts posterior draws of intercept and slope functions from saved draws of model parameters. A plot may be obtained if plot = TRUE. A list is returned invisibly with two fields.

beta.samp

a matrix with nmc many columns and (p+1)*length(tau.grid) many rows.

beta.est

a list of length (p+1), j-th element giving a 3-column matrix of median, 2.5th and 97.5th percentiles of the posterior distribution of {beta}[j]

See Also

qrjoint and summary.qrjoint for model fitting under qrjoint. Also see getBands for plotting credible bands for coefficients.

Examples

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
 
## Plasma data analysis

# recoding variables
data(plasma)
plasma$Sex <- as.factor(plasma$Sex)
plasma$SmokStat <- as.factor(plasma$SmokStat)
plasma$VitUse <- 3 - plasma$VitUse
plasma$VitUse <- as.factor(plasma$VitUse)

# creating predictors and response (beta carotene concentration in the plasma)
X <- model.matrix(BetaPlasma ~ Age + Sex + SmokStat + Quetelet + VitUse + Calories + 
        Fat + Fiber + Alcohol + Cholesterol + BetaDiet, data = plasma)[,-1]
Y <- plasma$BetaPlasma

# model fitting with 50 posterior samples from 100 iterations (thin = 2)
fit.qrj <- qrjoint(X, Y, 50, 2)

## Not run: 
betas <- coef(fit.qrj) ## no plots

## End(Not run)
betas <- coef(fit.qrj, plot = TRUE, col = "darkgreen") ## estimates are plotted