Fract.Poly {CorrMixed}R Documentation

Fit fractional polynomials

Description

Fits regression models with mm terms of the form XpX^{p}, where the exponents pp are selected from a small predefined set SS of both integer and non-integer values.

Usage

Fract.Poly(Covariate, Outcome, S=c(-2,-1,-0.5,0,0.5,1,2,3), Max.M=5, Dataset)

Arguments

Covariate

The covariate to be considered in the models.

Outcome

The outcome to be considered in the models.

S

The set SS from which each power pmp^{m} is selected. Default S={-2, -1, -0.5, 0, 0.5, 1, 2, 3}.

Max.M

The maximum order MM to be considered for the fractional polynomial. This value can be 55 at most. When M=5M=5, then fractional polynomials of order 11 to 55 are considered. Default Max.M=5.

Dataset

A data.frame that should consist of multiple lines per subject ('long' format).

Value

Results.M1

The results (powers and AIC values) of the fractional polynomials of order 11.

Results.M2

The results (powers and AIC values) of the fractional polynomials of order 22.

Results.M3

The results (powers and AIC values) of the fractional polynomials of order 33.

Results.M4

The results (powers and AIC values) of the fractional polynomials of order 44.

Results.M5

The results (powers and AIC values) of the fractional polynomials of order 55.

Author(s)

Wim Van der Elst, Geert Molenberghs, Ralf-Dieter Hilgers, & Nicole Heussen

References

Van der Elst, W., Molenberghs, G., Hilgers, R., & Heussen, N. (2015). Estimating the reliability of repeatedly measured endpoints based on linear mixed-effects models. A tutorial. Submitted.

Examples

# Open data
data(Example.Data)

# Fit fractional polynomials, mox. order = 3
FP <- Fract.Poly(Covariate = Time, Outcome = Outcome, 
Dataset = Example.Data, Max.M=3)

# Explore results
summary(FP)
# best fitting model (based on AIC) for m=3, 
# powers:  p_{1}=3,  p_{2}=3, and  p_{3}=2

# Fit model and compare with observed means

    # plot of mean 
Spaghetti.Plot(Dataset = Example.Data, Outcome = Outcome,
Time = Time, Id=Id, Add.Profiles = FALSE, Lwd.Me=1, 
ylab="Mean Outcome")

    # Coding of predictors (note that when p_{1}=p_{2},
    # beta_{1}*X ** {p_{1}} + beta_{2}*X ** {p_{1}} * log(X)
    #  and when p=0, X ** {0}= log(X)    )
term1 <- Example.Data$Time**3
term2 <- (Example.Data$Time**3) * log(Example.Data$Time) 
term3 <- Example.Data$Time**2 

    # fit model
Model <- lm(Outcome~term1+term2+term3, data=Example.Data)
Model$coef # regression weights (beta's)  

    # make prediction for time 1 to 47
term1 <- (1:47)**3 
term2 <- ((1:47)**3) * log(1:47)   
term3 <- (1:47)**2 
    # compute predicted values
pred <- Model$coef[1] + (Model$coef[2] * term1) + 
  (Model$coef[3] * term2) + (Model$coef[4] * term3)
   # Add predicted values to plot
lines(x = 1:47, y=pred,  lty=2)
legend("topright", c("Observed", "Predicted"), lty=c(1, 2))

[Package CorrMixed version 1.1 Index]