TABLE OF CONTENTS
methodsUni/unirec.formula [ Functions ]
NAME
unirec.formula --- method for formulas
FUNCTION
This is the main user-facing method of the package for univariate data. Its function is to conduct basic input checking and evaluate the formula, then construct a data frame in the format needed by bivrec.agdata for fitting. This entails creating a dummy second process, since bivrec.agdata works with bivariate processes. After bivrec.agdata is done, all references to the second dummy process are removed in cleanup.
SYNOPSIS
84 unirec.formula <- function(formula, data = parent.frame(), ...) 85 {
INPUTS
See package documentation
OUTPUTS
An object of class unirec. See package documentation
SOURCE
88 # in part based on coxph function 89 # Some of the code here is duplicated from bivrec.formula, but I don't see 90 # a good way around it 91 call <- match.call() 92 m <- match.call(expand.dots = FALSE) 93 if(is.matrix(eval(m$data, sys.parent()))) m$data <- as.data.frame(data) 94 m$...<-NULL 95 m[[1]] <- as.name("model.frame") 96 special <- c("id", "cluster", "strata") 97 Terms <- if (missing(data)) terms(formula, special) 98 else terms(formula, special, data = data) 99 m$formula <- Terms 100 m <- eval(m, sys.parent()) 101 n <- nrow(m) 102 resp <- model.extract(m, "response") 103 if (!is.Surv(resp)) stop("model response must be a Surv object") 104 if(dim(resp)[2] != 3) stop("response must be of type Surv(start, stop, status)") 105 start <- resp[, 1] 106 stop <- resp[, 2] 107 delta <- resp[, 3] 108 dropx <- NULL 109 clusterind <- attr(Terms, "specials")$cluster 110 111 # Cluster handling 112 clusternames <- NULL 113 if(length(clusterind) > 0){ 114 cluster <- m[, clusterind] 115 if(is.factor(cluster)) clusternames <- levels(cluster) 116 if(is.numeric(cluster)) clusternames <- as.character(unique(cluster)) 117 i <- as.numeric(as.factor(cluster)) 118 tempc <- untangle.specials(Terms, "cluster", 1:10) 119 ord <- attr(Terms, "order")[tempc$terms] 120 if (any(ord > 1)) stop("Cluster cannot be used in an interaction") 121 dropx <- c(dropx, tempc$terms) 122 }else{ 123 cluster <- rep(1, n) 124 i <- rep(1, n) 125 } 126 127 # ID handling 128 idind <- attr(Terms, "specials")$id 129 if(length(idind) == 0) stop("a subject id (unique within clusters) is required") 130 else{ 131 id <- m[, idind] 132 subjnames <- paste(cluster, id, sep = ".") 133 j <- rep(0, n) 134 for(ii in 1:max(i)){ 135 iids <- unique(id[i == ii]) 136 jj <- 1 137 for(thisid in iids){ 138 j[i == ii & id == thisid] <- jj 139 jj <- jj + 1 140 } 141 } 142 tempi <- untangle.specials(Terms, "id", 1:10) 143 ord <- attr(Terms, "order")[tempi$terms] 144 if (any(ord > 1)) stop("id cannot be used in an interaction") 145 dropx <- c(dropx, tempi$terms) 146 } 147 148 # Stratum handling 149 stratind <- attr(Terms, "specials")$strata 150 if(length(stratind) > 0) 151 { 152 strat <- m[, stratind] 153 ustrat <- unique(strat) 154 stratnames <- as.character(ustrat) 155 r <- rep(0, n) 156 for(rr in 1:length(ustrat)){ 157 r[strat == ustrat[rr]] <- rr 158 } 159 tempr <- untangle.specials(Terms, "strata", 1:10) 160 ord <- attr(Terms, "order")[tempr$terms] 161 if (any(ord > 1)) stop("strata cannot be used in an interaction") 162 dropx <- c(dropx, tempr$terms) 163 }else{ 164 r <- rep(1, n) 165 stratnames <- "1" 166 } 167 168 # Construct data frame to pass to unirec.agdata 169 Ki <- table(i * 1e6 + j) 170 k <- unlist(as.vector(sapply(Ki, function(x) 1:x))) 171 newTerms <- if(length(dropx)) Terms[ - dropx] else Terms 172 X <- model.matrix(newTerms, m) 173 X <- X[, -1, drop = FALSE] 174 agdata <- as.data.frame(cbind(i, j, k, r, start, stop, delta, X)) 175 sortord <- order(agdata$i, agdata$j) 176 agdata <- agdata[sortord, ] 177 subjnames <- unique(subjnames[sortord]) 178 processname <- as.character(call$formula[[2]][[4]]) 179 fit <- unirec(agdata, ...) 180 fit$call <- call 181 fit <- cleanunirecoutput(fit, clusternames, subjnames, stratnames, processname) 182 return(fit) 183 }