version 1.0.1

DESCRIPTION

@@ -0,0 +1,25 @@
+Package: OnAge
+Title: Test of Between-Group Differences in the Onset of Senescence
+Version: 1.0.1
+Date: 2017-10-20
+Author: Laurent Jacob, Frédéric Douhard, Jean-François Lemaître, Jean-Michel Gaillard, Aurélie Siberchicot
+Maintainer: Aurélie Siberchicot <aurelie.siberchicot@univ-lyon1.fr>
+Description: Implementation of a likelihood ratio test of differential
+    onset of senescence between two groups. Given two groups with
+    measures of age and of an individual trait likely to be subjected to
+    senescence (e.g. body mass), 'OnAge' provides an asymptotic p-value
+    for the null hypothesis that senescence starts at the same age in
+    both groups. The package implements the procedure used in 
+    Douhard et al. (2017) <doi:10.1111/oik.04421>.
+URL: https://lbbe.univ-lyon1.fr/OnAge.html
+License: GPL-3
+LazyLoad: yes
+Imports: stats
+Suggests: lme4
+Depends: R (>= 3.0.0)
+NeedsCompilation: no
+BuildVignettes: true
+Encoding: UTF-8
+Packaged: 2017-10-20 09:35:34 UTC; aurelie
+Repository: CRAN
+Date/Publication: 2017-10-20 10:21:02 UTC

MD5

@@ -0,0 +1,16 @@
+584b36f51a2e1bbef60f147b7f07002e *DESCRIPTION
+92c82e3e80fcde0af7ea5b4aaafa34fb *NAMESPACE
+3ea8d1a62f545a9828e659a5b65078b6 *NEWS
+3fabb10b6a196101882863f1c819ad71 *R/onset.test.R
+77a9146412af5eb47a193bd044d268d2 *build/vignette.rds
+dcc57255973abd1955435d46d3642e13 *data/RoeDeerMassData.RData
+61c205331835aef7cbeeea66587ff5d4 *inst/CITATION
+d32239bcb673463ab874e80d47fae504 *inst/COPYING
+afc0db49b6ffe87f76635acdcfd2cf78 *inst/doc/RoeDeerMass-simulation.R
+bc5ab8a7598ea5096d9501c515495fb5 *inst/doc/RoeDeerMass-simulation.Rnw
+ca51ab46ee1e75ac40ec77aa67f01320 *inst/doc/RoeDeerMass-simulation.pdf
+826c6977faf7a2f7d2900ce43563918c *man/RoeDeerMassData.Rd
+6c7a088cf16e46a70777aaac18a06945 *man/onset.test.Rd
+bc5ab8a7598ea5096d9501c515495fb5 *vignettes/RoeDeerMass-simulation.Rnw
+c34c2b384e0887ededcf4cba46357747 *vignettes/bibli.bib
+59edcbef45afebb72c9739f8abf2c529 *vignettes/pre-computed-simulation.RData

NAMESPACE

@@ -0,0 +1,5 @@
+export('onset.test')
+
+importFrom('grDevices', 'dev.off', 'pdf')
+importFrom('graphics', 'abline', 'par', 'plot', 'points')
+importFrom('stats', 'optimize', 'pchisq', 'qchisq', 'uniroot')

NEWS

@@ -0,0 +1,7 @@
+* OnAge 1.0-1: First release
+
+* August 28th 2017: created. 'OnAge' implements a likelihood ratio test
+of differential onset on senescence. Given two populations with
+measured ages and a senescence factor (e.g. body mass), 'OnAge' provides
+an asymptotic p-value for the null hypothesis that senescence starts
+at the same age in both populations.

R/onset.test.R

@@ -0,0 +1,133 @@
+## Build a likelihood ratio test of H0 (same onset in both groups) vs
+## H1 (different onset allowed)
+onset.test <- function(ll, data1, data2, search.range, CI.lvl=0.95, tol=.Machine$double.eps^0.25, warn=FALSE, do.plot=FALSE, plot.file=NULL, grid.len=100){
+    ## Maximize the likelihood over all parameters including onset for
+    ## group 1
+    c1.opt <- optimize(ll, interval=search.range, dataIn = data1, maximum=TRUE, tol=tol)
+
+    ## Maximize the likelihood over all parameters including onset for
+    ## group 2
+    c2.opt <- optimize(ll, interval=search.range, dataIn = data2, maximum=TRUE, tol=tol)
+
+    ## Joint loglikelihood across data1 and data2 with a single onset (ie,
+    ## loglikelihood under H0)
+    ll.joint <- function(thr, ll, data1, data2){
+        ll(thr, data1) + ll(thr, data2)
+    }
+    ## Maximize the likelihood over all parameters including (single)
+    ## onset for both groups    
+    joint.opt <- optimize(ll.joint, interval=search.range, ll, data1 = data1, data2 = data2, maximum=TRUE, tol=tol)
+
+    ## Loglikelihood under H1 (sum of the loglikelihood obtained
+    ## for both groups)    
+    lh1 <- as.numeric(c1.opt$objective + c2.opt$objective)
+    ## Loglikelihood under H0
+    lh0 <- as.numeric(joint.opt$objective)
+
+    ## Loglikelihood ratio (LLR) statistic
+    llr <- 2*(lh1 - lh0)
+
+    ## Sanity check: the likelihood under H1 should be larger than
+    ## under H0 by construction.    
+    cvg.ok <- TRUE
+    if(llr < 0){
+        ## lh1 <- lh0
+        llr <- 0
+        if(warn){
+            warning('[onset.test] Likelihood under H1 was lower than under H0, probably due to optimization issues. Returning null log likelihood ratio. est.1 and est.2 values are suboptimal and should be both replaced by est.joint.')
+        }
+        cvg.ok <- FALSE
+    }
+
+    ## LLR is asymptotically chi2_1 under H0. Use this to compute an
+    ## (asymptotic) p-value.    
+    pv <- 1 - pchisq(llr, 1)
+
+    ## Optionally get confidence intervals on the ages at onset
+    if(!is.na(CI.lvl)){
+      CI.1 <- CI.from.opt(c1.opt, ll, lvl=CI.lvl, search.range, data1)
+      CI.2 <- CI.from.opt(c2.opt, ll, lvl=CI.lvl, search.range, data2)
+      joint.CI <- CI.from.opt(joint.opt, ll.joint, lvl=CI.lvl, search.range, ll, data1, data2)
+    }else{
+        CI.1 <- CI.2 <- joint.CI <- NA
+    }
+
+    if(do.plot){
+
+        ll.plot <- function(ll.fun, ..., x.grid, x.opt, ll.opt, CI, data.name, CEX=1.5){
+            ll.grid <- sapply(x.grid, ll.fun, ...)
+            plot(x.grid, ll.grid, main=data.name, xlab="Age at onset", ylab="Log likelihood", cex.lab=CEX, cex.axis=CEX)
+            points(x.opt, ll.opt, pch=8, col='red')
+            if(!any(is.na(CI))){
+                abline(v=CI[1], lty=2)
+                abline(v=CI[2], lty=2)
+                abline(h=min(ll.fun(CI[1], ...), ll.fun(CI[2], ...)), lty=2)
+            }
+        }
+
+        if(!is.null(plot.file)){
+            pdf(file=plot.file, width=15, height=5)
+        }
+        mars <- c(c(4.5, 5, 2.5, 1))
+        par(mar=mars, mfrow=c(1, 3))
+        thr.grid <- seq(search.range[1], search.range[2], length.out=grid.len)
+        ll.plot(ll, data1, x.grid=thr.grid, x.opt=c1.opt$maximum, ll.opt=c1.opt$objective, CI=CI.1, data.name='Group 1')
+        ll.plot(ll, data2, x.grid=thr.grid, x.opt=c2.opt$maximum, ll.opt=c2.opt$objective, CI=CI.2, data.name='Group 2')
+        ll.plot(ll.joint, ll, data1, data2, x.grid=thr.grid, x.opt=joint.opt$maximum, ll.opt=joint.opt$objective, CI=joint.CI, data.name='All data')
+
+        if(!is.null(plot.file)){
+            dev.off()
+        }        
+    }
+
+    return(list(pv=pv,
+                est.1=c1.opt$maximum, est.2=c2.opt$maximum, est.joint=joint.opt$maximum,
+                CI.1=CI.1, CI.2=CI.2, joint.CI=joint.CI,
+                lh0=lh0, lh1=lh1, llr=llr,
+                cvg.ok=cvg.ok))
+}
+
+
+## Build CI by identiying smallest and largest values thr of onset for
+## which the null hypothesis true_onset=thr is rejected at level
+## 1-lvl. Caveat: assumes the log-likelihood is monotonically
+## decreasing, i.e., that the log-likelihood ratio statistic is
+## monotonically increasing around the maximum likelihood estimator of
+## onset (can be checked visually using do.plot=TRUE in
+## onset.test). Otherwise the CI may not be connected and a more
+## sophisticated method should be used.
+
+CI.from.opt <- function(opt, ll.fun, lvl=0.95, search.range, ...){
+
+    ## Zero of this function is an onset thresholds thr for which
+    ## H0:true_onset=thr is rejected at level exactly 1-lvl using a
+    ## likelihood ratio test. The ML estimator of onset leads to
+    ## rejection at level 1 by construction of the likelihood ratio
+    ## statistic (equal to 0).  We assume (monotonicity of the
+    ## likelihood) that all values between the ML estimator and this
+    ## zero lead to rejection at level >= 1-lvl, making them part of
+    ## the lvl confidence interval.    
+    dist.from.quantile <- function(thr)
+    {
+        2*(opt$objective - ll.fun(thr, ...)) - qchisq(lvl, 1)
+    }
+
+    bp.opt <- opt$maximum # Maximum likelihood onset
+
+    if(dist.from.quantile(search.range[1]) < 0){
+        ## warning(sprintf('[CI.from.opt] H0: onset=search.range[1] not rejected at level %g. Setting lower bound of the CI to search.range[1].', 1-lvl))
+        lb <- search.range[1]
+    }else{
+        lb <- uniroot(dist.from.quantile, lower=search.range[1], upper=bp.opt)$root
+    }
+    if(dist.from.quantile(search.range[2]) < 0){
+        ## warning(sprintf('[CI.from.opt] H0: onset=search.range[2] not rejected at level %g. Setting upper bound of the CI to search.range[2].', 1-lvl))
+        ub <- search.range[2]
+    }else{
+        ub <- uniroot(dist.from.quantile, lower=bp.opt, upper=search.range[2])$root
+    }
+    CI <- c(lb, ub)    
+    names(CI) <- c('lb', 'ub')
+
+    return(CI)
+}

inst/CITATION

@@ -0,0 +1,25 @@
+citEntry(
+  # - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 
+  # BibTeX entry:
+  # - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
+    entry = "article",
+    author = "Douhard, Frederic and Gaillard, Jean-Michel and Pellerin, Maryline and Jacob, Laurent and Lemaitre, Jean-Francois",
+    title = "The cost of growing large: costs of post-weaning growth on body mass senescence in a wild mammal",
+    journal = "Oikos",
+    volume = "126",
+    number = "9",
+    pages = "1329--1338",
+    publisher = "Blackwell Publishing Ltd",
+    issn = "1600-0706",
+    url = "http://dx.doi.org/10.1111/oik.04421",
+    doi = "10.1111/oik.04421",
+    year = "2017",
+  # - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 
+  # Plain-text citation:
+  # - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 
+    textVersion = paste("Douhard, F., Gaillard, J.-M., Pellerin, M., Jacob, L. and Lemaître, J.-F. (2017),",
+                        "The cost of growing large: costs of post-weaning growth on body mass senescence in a wild mammal.",
+                        "Oikos, 126: 1329–1338.", 
+                        "doi:10.1111/oik.04421"
+                        )
+)