# Toy Life History Data having Directions of Recession Toy life history data created to exhibit the phenomenon of directions of recession. It was analyzed in a special topics course on aster models Fall 2018. data ```r data(foobar) ``` ## Format `data(foobar)` loads four R objects. - **fam**: a vector of family indices. - **pred**: a vector of predecessor indices. - **vars**: a vector of variable names associated with the nodes of the aster graph. - **redata**: a data frame having data on 300 individuals, each of which has `length(vars) == 4` components of fitness, so the aster graph for one individual has four nodes. This data frame is already in long format; no need to `reshape`. The variables in this data frame are: - **resp**: the response vector. - **varb**: Categorical. Gives node of graphical model corresponding to each component of `resp`. See details below. - **root**: All ones. Root variables for graphical model. - **id**: Categorical. Individual ID numbers. - **trt**: Categorical. Treatment. - **blk**: Categorical. Block. - **fit**: Bernoulli (zero-or-one-valued). Indicator variable of the fitness nodes of the graph; in these data there is just one node for fitness. ## Details The levels of `varb` indicate nodes of the graphical model to which the corresponding elements of the response vector `resp` belong. This is the typical long format produced by the R `reshape` function. For each individual, there are several response variables. All response variables are combined in one vector `resp`. The variable `varb` indicates which original variable the the corresponding component of the response vector was. The variable `id` indicates which individual the corresponding component of the response vector was. ## Source Charles J. Geyer [https://www.stat.umn.edu/geyer/8931aster/foobar.rda](https://www.stat.umn.edu/geyer/8931aster/foobar.rda) ## References These data are analyzed in deck 9 of the slides for a special topics course on aster models taught fall semester 2018 ([https://www.stat.umn.edu/geyer/8931aster/slides/s9.pdf](https://www.stat.umn.edu/geyer/8931aster/slides/s9.pdf)). ## Examples ```r data(foobar) library(aster) aout <- aster(resp ~ varb + fit : (trt * blk), pred, fam, varb, id, root, data = redata) foo <- try(summary(aout)) # gives an error about directions of recession # get directions of recession dor <- attr(foo, "condition")$dor dor # found one apparent direction of recession # from regular pattern # it looks like a true direction of recession dor <- dor / max(abs(dor)) dor # but what does it do? For that need to map to saturated model # parameter space modmat <- aout$modmat dim(modmat) # oof! modmat is three-dimensional. Need an actual matrix here. modmat <- matrix(as.vector(modmat), ncol = length(dor)) dor.phi <- drop(modmat %*% dor) names(dor.phi) <- with(redata, paste(id, as.character(varb), sep = ".")) dor.phi[dor.phi != 0] fam.default()[fam[vars == "seeds"]] # since the support of the Poisson distribution is bounded above, # actually this must be minus the DOR (if it is a DOR at all). # check that all components of response vector for which dor.phi == 1 are zero # (lower bound of Poisson range) all(redata$resp[dor.phi == 1] == 0) # so minus dor.phi is a true direction of recession in the saturated model # canonical parameter space, and minus dor is a true direction of recession # in the submodel canonical parameter space # # try to get more info trt.blk <- with(redata, paste(as.character(trt), as.character(blk), sep = ".")) unique(trt.blk[dor.phi == 1]) # the reason for the direction of recession is that every individual getting # treatment a in block A had zero seeds. # # the reason the submodel DOR, R object dor, was so hard to interpret was # because fit:trta:blkA is not in the model. So let's force it in redata <- transform(redata, trt = relevel(trt, ref = "b"), blk = relevel(blk, ref = "B")) # Note: following code is copied exactly from above. Only difference # is releveling in the immediately preceding statement aout <- aster(resp ~ varb + fit : (trt * blk), pred, fam, varb, id, root, data = redata) foo <- try(summary(aout)) dor <- attr(foo, "condition")$dor dor <- dor / max(abs(dor)) dor # now it is obvious from looking at this dor that all individuals in trt a # and blk A are at the lower end (zero) of the Poisson range. # maybe the other dor we had previously would be "obvious" to someone # sufficiently skilled in understanding the meaning of regression coefficients # but not "obvious" to your humble author # # as for what to do about this, see the course slides cited in the reference # section. There is no single Right Thing to do. ```

foobar dataset