diff -Nru r-cran-statmod-1.4.30/debian/changelog r-cran-statmod-1.4.32/debian/changelog --- r-cran-statmod-1.4.30/debian/changelog 2018-06-22 21:55:41.000000000 +0000 +++ r-cran-statmod-1.4.32/debian/changelog 2019-07-02 12:16:07.000000000 +0000 @@ -1,3 +1,14 @@ +r-cran-statmod (1.4.32-1) unstable; urgency=medium + + * New upstream version 1.4.32 + * Update long package description + * Bump to debhelper compat level 12 + * Add Testsuite: autopkgtest-pkg-r + * Add Rules-Requires-Root: no + * Bump to S-V 4.3.0 + + -- Sébastien Villemot Tue, 02 Jul 2019 14:16:07 +0200 + r-cran-statmod (1.4.30-2) unstable; urgency=medium * Team upload. diff -Nru r-cran-statmod-1.4.30/debian/compat r-cran-statmod-1.4.32/debian/compat --- r-cran-statmod-1.4.30/debian/compat 2018-06-22 21:55:41.000000000 +0000 +++ r-cran-statmod-1.4.32/debian/compat 1970-01-01 00:00:00.000000000 +0000 @@ -1 +0,0 @@ -11 diff -Nru r-cran-statmod-1.4.30/debian/control r-cran-statmod-1.4.32/debian/control --- r-cran-statmod-1.4.30/debian/control 2018-06-22 21:55:41.000000000 +0000 +++ r-cran-statmod-1.4.32/debian/control 2019-07-02 12:15:54.000000000 +0000 @@ -3,13 +3,15 @@ Uploaders: Sébastien Villemot Section: gnu-r Priority: optional -Build-Depends: debhelper (>= 11~), +Build-Depends: debhelper-compat (= 12), dh-r, r-base-dev -Standards-Version: 4.1.4 +Standards-Version: 4.3.0 Vcs-Browser: https://salsa.debian.org/r-pkg-team/r-cran-statmod Vcs-Git: https://salsa.debian.org/r-pkg-team/r-cran-statmod.git Homepage: https://cran.r-project.org/package=statmod +Testsuite: autopkgtest-pkg-r +Rules-Requires-Root: no Package: r-cran-statmod Architecture: any @@ -20,9 +22,9 @@ Suggests: ${R:Suggests} Description: GNU R package providing algorithms and functions for statistical modeling This R package provides a collection of algorithms and functions to aid - statistical modeling. It includes growth curve comparisons, limiting dilution - analysis (aka ELDA), mixed linear models, heteroscedastic regression, + statistical modeling. Includes limiting dilution analysis (aka ELDA), growth + curve comparisons, mixed linear models, heteroscedastic regression, inverse-Gaussian probability calculations, Gauss quadrature and a secure - convergence algorithm for nonlinear models. It also includes advanced - generalized linear model functions that implement secure convergence, - dispersion modeling and Tweedie power-law families. + convergence algorithm for nonlinear models. Also includes a number of advanced + generalized linear model functions including new Tweedie and Digamma glm + families and a secure convergence algorithm. \ No newline at end of file diff -Nru r-cran-statmod-1.4.30/DESCRIPTION r-cran-statmod-1.4.32/DESCRIPTION --- r-cran-statmod-1.4.30/DESCRIPTION 2017-06-18 14:56:34.000000000 +0000 +++ r-cran-statmod-1.4.32/DESCRIPTION 2019-05-29 09:40:02.000000000 +0000 @@ -1,15 +1,15 @@ Package: statmod -Version: 1.4.30 -Date: 2017-06-16 +Version: 1.4.32 +Date: 2019-05-29 Title: Statistical Modeling Author: Gordon Smyth [cre, aut], Yifang Hu [ctb], Peter Dunn [ctb], Belinda Phipson [ctb], Yunshun Chen [ctb] Maintainer: Gordon Smyth -Depends: R (>= 1.6.1) +Depends: R (>= 3.0.0) Imports: stats, graphics Suggests: MASS, tweedie -Description: A collection of algorithms and functions to aid statistical modeling. Includes growth curve comparisons, limiting dilution analysis (aka ELDA), mixed linear models, heteroscedastic regression, inverse-Gaussian probability calculations, Gauss quadrature and a secure convergence algorithm for nonlinear models. Includes advanced generalized linear model functions that implement secure convergence, dispersion modeling and Tweedie power-law families. +Description: A collection of algorithms and functions to aid statistical modeling. Includes limiting dilution analysis (aka ELDA), growth curve comparisons, mixed linear models, heteroscedastic regression, inverse-Gaussian probability calculations, Gauss quadrature and a secure convergence algorithm for nonlinear models. Also includes a number of advanced generalized linear model functions including new Tweedie and Digamma glm families and a secure convergence algorithm. License: GPL-2 | GPL-3 NeedsCompilation: yes -Packaged: 2017-06-17 23:53:15 UTC; smyth +Packaged: 2019-05-29 08:16:02 UTC; smyth Repository: CRAN -Date/Publication: 2017-06-18 14:56:34 UTC +Date/Publication: 2019-05-29 09:40:02 UTC diff -Nru r-cran-statmod-1.4.30/inst/NEWS r-cran-statmod-1.4.32/inst/NEWS --- r-cran-statmod-1.4.30/inst/NEWS 2017-06-16 09:03:18.000000000 +0000 +++ r-cran-statmod-1.4.32/inst/NEWS 2019-05-29 08:15:57.000000000 +0000 @@ -1,3 +1,15 @@ +29 May 2019: statmod 1.4.32 + +- Bug fix to glmnb.fit() when all the y values are zero. + +24 Feb 2017: statmod 1.4.31 + +- Add Dunn & Smyth (2018) GLM book as a reference for glm.scoretest() + and Randomized Quantile Residuals. Rewrite the Details section of + the glm.scoretest help page. + +- Update minimum version of R to 3.0.0 in DESCRIPTION file. + 16 June 2017: statmod 1.4.30 - Bug fix to qinvgauss(). In some case the gamma approximation used diff -Nru r-cran-statmod-1.4.30/man/fitNBP.Rd r-cran-statmod-1.4.32/man/fitNBP.Rd --- r-cran-statmod-1.4.30/man/fitNBP.Rd 2009-07-13 04:04:00.000000000 +0000 +++ r-cran-statmod-1.4.32/man/fitNBP.Rd 2019-05-29 08:13:12.000000000 +0000 @@ -34,6 +34,11 @@ This iteration is monotonically convergent for the dispersion. } +\note{ +This function has been made obsolete by the edgeR package on Bioconductor, see +\url{https://bioconductor.org/packages/edgeR}. +} + \value{ List with components \item{coefficients}{numeric matrix of rates for each tag (gene) and each group} @@ -53,6 +58,9 @@ \seealso{ \code{\link[statmod]{sage.test}} + +The edgeR package on Biconductor provides new and better functions to fit negative-binomial glms to SAGE or RNA-seq data. +See particularly the \code{glmFit} and \code{mglmOneWay} functions. } \examples{ diff -Nru r-cran-statmod-1.4.30/man/glmgam.Rd r-cran-statmod-1.4.32/man/glmgam.Rd --- r-cran-statmod-1.4.30/man/glmgam.Rd 2013-11-19 23:30:48.000000000 +0000 +++ r-cran-statmod-1.4.32/man/glmgam.Rd 2019-05-29 07:47:31.000000000 +0000 @@ -1,24 +1,16 @@ \name{glmgam.fit} \alias{glmgam.fit} -\alias{glmnb.fit} -\title{Fit Generalized Linear Model by Fisher Scoring with Levenberg Damping} +\title{Fit Gamma Generalized Linear Model by Fisher Scoring with Identity Link} \description{ Fit a generalized linear model with secure convergence. -Provided for gamma glm with identity links or negative binomial glm with log-links. } \usage{ -glmgam.fit(X, y, coef.start=NULL, tol=1e-6, maxit=50, trace=FALSE) -glmnb.fit(X, y, dispersion, weights=NULL, offset=0, coef.start=NULL, start.method="mean", - tol=1e-6, maxit=50, trace=FALSE) +glmgam.fit(X, y, coef.start = NULL, tol = 1e-6, maxit = 50, trace = FALSE) } \arguments{ \item{X}{design matrix, assumed to be of full column rank. Missing values not allowed.} \item{y}{numeric vector of responses. Negative or missing values not allowed.} - \item{dispersion}{numeric vector of over-dispersion parameters for negative binomial. If of length 1, then same over-dispersion is assumed for all observations.} - \item{weights}{numeric vector of positive weights, defaults to all one.} - \item{offset}{offset vector for linear model} \item{coef.start}{numeric vector of starting values for the regression coefficients} - \item{start.method}{method used to find starting values, possible values are \code{"mean"} and \code{"log(y)"}} \item{tol}{small positive numeric value giving convergence tolerance} \item{maxit}{maximum number of iterations allowed} \item{trace}{logical value. If \code{TRUE} then output diagnostic information at each iteration.} @@ -30,21 +22,32 @@ \item{deviance}{residual deviance} \item{iter}{number of iterations used to convergence. If convergence was not achieved then \code{iter} is set to \code{maxit+1}.} } + \details{ -These functions implement a modified Fisher scoring algorithm for generalized linear models, similar to the Levenberg-Marquardt algorithm for nonlinear least squares. +This function implements a modified Fisher scoring algorithm for generalized linear models, similar to the Levenberg-Marquardt algorithm for nonlinear least squares. The Levenberg-Marquardt modification checks for a reduction in the deviance at each step, and avoids the possibility of divergence. The result is a very secure algorithm that converges for almost all datasets. -\code{glmgam.fit} is in principle similar to \code{glm.fit(X,y,family=Gamma(link="identity"))} but with much more secure convergence. -This function is used by \code{\link{mixedModel2Fit}}. +\code{glmgam.fit} is in principle equivalent to \code{glm.fit(X,y,family=Gamma(link="identity"))} but with much more secure convergence. +} -\code{glmnb.fit} is in principle similar to \code{glm.fit(X,y,family=negative.binomial(link="log",theta=1/dispersion))} but with more secure convergence. +\author{Gordon Smyth and Yunshun Chen} + +\references{ +Dunn, PK, and Smyth, GK (2018). +\emph{Generalized linear models with examples in R}. Springer, New York, NY. +\url{https://doi.org/10.1007/978-1-4419-0118-7} +} +\seealso{ +\code{glmgam.fit} is called by \code{\link{mixedModel2Fit}}. + +\code{\link{glm}} is the standard glm fitting function in the stats package. } -\author{Gordon Smyth and Yunshun Chen} + \examples{ -y <- rgamma(10,shape=5) -X <- cbind(1,1:10) -fit <- glmgam.fit(X,y,trace=TRUE) +y <- rgamma(10, shape=5) +X <- cbind(1, 1:10) +fit <- glmgam.fit(X, y, trace=TRUE) } \keyword{regression} diff -Nru r-cran-statmod-1.4.30/man/glmnbfit.Rd r-cran-statmod-1.4.32/man/glmnbfit.Rd --- r-cran-statmod-1.4.30/man/glmnbfit.Rd 1970-01-01 00:00:00.000000000 +0000 +++ r-cran-statmod-1.4.32/man/glmnbfit.Rd 2019-05-29 08:04:41.000000000 +0000 @@ -0,0 +1,60 @@ +\name{glmnb.fit} +\alias{glmnb.fit} +\title{Fit Negative Binomial Generalized Linear Model with Log-Link} +\description{ +Fit a generalized linear model with secure convergence. +} +\usage{ +glmnb.fit(X, y, dispersion, weights = NULL, offset = 0, coef.start = NULL, + start.method = "mean", tol = 1e-6, maxit = 50L, trace = FALSE) +} +\arguments{ + \item{X}{design matrix, assumed to be of full column rank. Missing values not allowed.} + \item{y}{numeric vector of responses. Negative or missing values not allowed.} + \item{dispersion}{numeric vector of dispersion parameters for the negative binomial distribution. If of length 1, then the same dispersion is assumed for all observations.} + \item{weights}{numeric vector of positive weights, defaults to all one.} + \item{offset}{offset vector for linear model} + \item{coef.start}{numeric vector of starting values for the regression coefficients} + \item{start.method}{method used to find starting values, possible values are \code{"mean"} or \code{"log(y)"}} + \item{tol}{small positive numeric value giving convergence tolerance} + \item{maxit}{maximum number of iterations allowed} + \item{trace}{logical value. If \code{TRUE} then output diagnostic information at each iteration.} +} +\value{ +List with the following components: + \item{coefficients}{numeric vector of regression coefficients} + \item{fitted}{numeric vector of fitted values} + \item{deviance}{residual deviance} + \item{iter}{number of iterations used to convergence. If convergence was not achieved then \code{iter} is set to \code{maxit+1}.} +} +\details{ +This function implements a modified Fisher scoring algorithm for generalized linear models, analogous to the Levenberg-Marquardt algorithm for nonlinear least squares. +The Levenberg-Marquardt modification checks for a reduction in the deviance at each step, and avoids the possibility of divergence. +The result is a very secure algorithm that converges for almost all datasets. + +\code{glmnb.fit} is in principle equivalent to \code{glm.fit(X,y,family=negative.binomial(link="log",theta=1/dispersion))} but with more secure convergence. +Here \code{negative.binomial} is a function in the MASS package. + +The \code{dispersion} parameter is the same as \code{1/theta} for the \code{MASS::negative.binomial} function or \code{1/size} for the \code{stats::rnbinom} function. +\code{dispersion=0} corresponds to the Poisson distribution. +} +\author{Gordon Smyth and Yunshun Chen} + +\references{ +Dunn, PK, and Smyth, GK (2018). +\emph{Generalized linear models with examples in R}. Springer, New York, NY. +\url{https://doi.org/10.1007/978-1-4419-0118-7} +} + +\seealso{ +The \code{glmFit} function in the edgeR package on Bioconductor is a high-performance version of \code{glmnb.fit} for many \code{y} vectors at once. + +\code{\link{glm}} is the standard glm fitting function in the stats package. +\code{negative.binomial} in the MASS package defines a negative binomial family for use with \code{glm}. +} +\examples{ +y <- rnbinom(10, mu=1:10, size=5) +X <- cbind(1, 1:10) +fit <- glmnb.fit(X, y, dispersion=0.2, trace=TRUE) +} +\keyword{regression} diff -Nru r-cran-statmod-1.4.30/man/glmscoretest.Rd r-cran-statmod-1.4.32/man/glmscoretest.Rd --- r-cran-statmod-1.4.30/man/glmscoretest.Rd 2010-01-05 07:01:30.000000000 +0000 +++ r-cran-statmod-1.4.32/man/glmscoretest.Rd 2019-05-29 07:33:50.000000000 +0000 @@ -17,27 +17,49 @@ } \details{ -Rao's score statistic. -Is the locally most powerful test for testing vs a one-sided alternative. -Asympotically equivalent to likelihood ratio tests, but convenient for one-sided tests. - -This function computes a score test statistics for adding each covariate individually. +Rao's score test is a type of asymptotic test that is an alternative to Wald tests or likelihood ratio tests (LRTs) (Dunn and Smyth, 2018). +Wald tests are computed by dividing parameter estimates by their standard errors. +LRTs are computed from differences in the log-likihoods between the null and alternative hypotheses. +Score tests are computed from log-likelihood derivatives. +All three types of tests (Wald, score and LRT) are asymptotically equivalent under ideal circumstances, but the score and LRT tests are invariant under-reparametrization whereas Wald tests are not. + +One of the main differences between the tests is the need for estimation of parameters under the null and alternative hypotheses. +Wald tests require maximum likelihood estimates (MLEs) to be computed only under the alternative hypothesis, LRT tests require both the null and alternative models to be fitted, while score tests require only the null hypothesis to be fitted. + +When a generalized linear model (GLM) is fitted in R using the \code{glm} function, the \code{summary()} function presents Wald tests for all the coefficients in the linear model while \code{anova()} is able to compute likelihood ratio tests. +GLM output in R has historically not included score tests, although score tests can be a very computationally coefficient choice when one wants to test for many potential additional covariates being added to a relatively simple null model. + +A number of well known Pearson chisquare statistics, including goodness of fit statistics and the Pearson test for independence in a contingency table can be derived as score tests (Smyth, 2003; Dunn and Smyth, 2018). + +This function computes score test statistics for adding a single numerical covariate to a GLM, given the \code{glm} output for the null model. +It makes very efficient use of the quantities already stored in the GLM fit object. +A computational formula for the score test statistics is given in Section 7.2.6 of Dunn and Smyth (2018). The dispersion parameter is treated as for \code{\link{summary.glm}}. -If \code{NULL}, the Pearson estimator is used, except for the binomial and Poisson +If \code{NULL}, the Pearson estimator is used, except for the binomial, Poisson and negative binomial families, for which the dispersion is one. + +Note that the \code{anova.glm} function in the stats package has offered a Rao score test option since 2011, but it requires model fits under the alternative as well as the null hypothesis, which does not take advantage of the computational advantage of score test. +The \code{glm.scoretest} is more efficient as it does not require a full model fit. +On the other hand, \code{anova.glm} can compute score tests for factors and multiple covariates, which \code{glm.scoretest} does not currently do. } -\value{numeric vector containing the z-statistics, one for each covariate.} +\value{ +numeric vector containing the z-statistics, one for each covariate. +The z-statistics can be treated as standard normal under the null hypothesis. +} \author{Gordon Smyth} \seealso{ -\code{\link{glm}}, \code{\link{add1}} +\code{\link{glm}}, \code{\link{add1}}, \code{\link{anova.glm}} } \references{ -Lovison, G (2005). On Rao score and Pearson $X^2$ statistics in generalized linear models. +Dunn, PK, and Smyth, GK (2018). \emph{Generalized linear models with examples in R}. Springer, New York, NY. +\url{https://doi.org/10.1007/978-1-4419-0118-7} + +Lovison, G (2005). On Rao score and Pearson X^2 statistics in generalized linear models. \emph{Statistical Papers}, 46, 555-574. Pregibon, D (1982). Score tests in GLIM with applications. @@ -55,15 +77,15 @@ # First the usual test y <- c(20,40,40,30) -chisq.test(matrix(y,2,2),correct=FALSE) +chisq.test(matrix(y,2,2), correct=FALSE) # Now same test using glm.scoretest a <- gl(2,1,4) b <- gl(2,2,4) -fit <- glm(y~a+b,family=poisson) +fit <- glm(y~a+b, family=poisson) x2 <- c(0,0,0,1) -z <- glm.scoretest(fit,x2) +z <- glm.scoretest(fit, x2) z^2 } diff -Nru r-cran-statmod-1.4.30/man/qresiduals.Rd r-cran-statmod-1.4.32/man/qresiduals.Rd --- r-cran-statmod-1.4.30/man/qresiduals.Rd 2009-09-29 02:41:41.000000000 +0000 +++ r-cran-statmod-1.4.32/man/qresiduals.Rd 2019-05-29 07:33:25.000000000 +0000 @@ -56,6 +56,10 @@ \references{ Dunn, K. P., and Smyth, G. K. (1996). Randomized quantile residuals. \emph{Journal of Computational and Graphical Statistics} \bold{5}, 1-10. \url{http://www.statsci.org/smyth/pubs/residual.html} + +Dunn, PK, and Smyth, GK (2018). +\emph{Generalized linear models with examples in R}. Springer, New York, NY. +\url{https://doi.org/10.1007/978-1-4419-0118-7} } \author{Gordon Smyth} diff -Nru r-cran-statmod-1.4.30/man/sage.test.Rd r-cran-statmod-1.4.32/man/sage.test.Rd --- r-cran-statmod-1.4.30/man/sage.test.Rd 2013-09-27 08:36:48.000000000 +0000 +++ r-cran-statmod-1.4.32/man/sage.test.Rd 2019-05-29 08:11:33.000000000 +0000 @@ -3,11 +3,6 @@ \title{Exact Binomial Tests For Comparing Two SAGE Libraries (Obsolete)} \description{ -This function is kept here so as not to break code that depends on it, -but has been replaced by \code{binomTest} in the edgeR Bioconductor package and -is no longer updated. -It may be removed in a later release of this package. - Computes p-values for differential abundance for each tag between two digital libraries, conditioning on the total count for each tag. The counts in each group as a proportion of the whole are assumed to follow a binomial distribution. @@ -47,6 +42,13 @@ This function has been made obsolete by \code{binomTest} in the edgeR package. } +\note{ +This function is kept in the statmod package so as not to break code that depends on it +but it has been replaced by \code{binomTest} in the edgeR Bioconductor package and +is no longer updated. +It may be removed in a later release of this package. +} + \value{ Numeric vector of p-values. } @@ -64,7 +66,10 @@ \author{Gordon Smyth} \seealso{ -\code{binomTest} (edgeR package), \code{\link{binom.test}} (stats package) +The \code{binomTest} function in the edgeR package on Bioconductor is a newer and better version of this function, see +\url{https://bioconductor.org/packages/edgeR}. + +\code{\link{binom.test}} in the stats package performs univariate binomial tests. } \examples{ diff -Nru r-cran-statmod-1.4.30/MD5 r-cran-statmod-1.4.32/MD5 --- r-cran-statmod-1.4.30/MD5 2017-06-18 14:56:34.000000000 +0000 +++ r-cran-statmod-1.4.32/MD5 2019-05-29 09:40:02.000000000 +0000 @@ -1,13 +1,13 @@ -2f4dabd9e0b70eedfd4ef8a82fda44cb *DESCRIPTION -e1f88fb7fb33316042315e919b39cd7d *NAMESPACE +6b65c94fe83cf0e8e87beeccd6739046 *DESCRIPTION +f9339c225b76e77744e0ff59cf229454 *NAMESPACE cabcedc80916e132755147508dc1808b *R/digamma.R b67d4c729d215c0952dc4594882da6fa *R/digammaf.R 2ca5d21ceedcedb7a85e9f46800187c7 *R/elda.R f48cdb780b1cc86bfeef8c4722912c2f *R/fitNBP.R 27838c3de5b6c7ff37ea6e92870903cc *R/forward.R 45cf210a8910c61246405ba70a0f1ae8 *R/gaussquad.R -7871d84512b208401cec51e401c990d9 *R/glmgam.R -99187037b66bc90438229911a23f1f2d *R/glmnb.R +7e07d1b75eaa61131248c0fd176ed84a *R/glmgam.R +e0a5de6d8c194cdd1d59907eca7b8370 *R/glmnb.R 87e2551f59f03a67cf4b0a9bf6479972 *R/glmscoretest.R 98e67ef55b5ac6c79e18608e6e6e2b91 *R/growthcurve.R 1747893341300d18e4af7032e4f91fc0 *R/hommel.R @@ -24,15 +24,16 @@ afad3149539912e729cb6d0aea8e53df *R/tweedief.R e600a474fedf570b1913b278022d46bf *data/welding.rdata 12b53c3de6776dbc33ac3a0e098f033a *inst/CITATION -997902dc06f56c32cc0b34441e4f2334 *inst/NEWS +d60c6df657a2b74606ad0533b35d0f16 *inst/NEWS d37cd00f740faf6738e9a0c992f54eda *man/digammaf.Rd 1ee58cb70d97c57b1aeb31e2f94258a9 *man/elda.Rd -2f065b8fef25cc88b80bd996eb24b0f1 *man/fitNBP.Rd +c2270a376b7575b1977f6d2587d00f56 *man/fitNBP.Rd 3508ac5839f19813b7894387b7a21bb7 *man/forward.Rd 9e4d97cd65ba17014c30878b39a9a351 *man/gauss.quad.Rd 6fdf80d724c6e363280eecb45419977d *man/gauss.quad.prob.Rd -809ccf2ff24e804a2a5bb17baa1c80b9 *man/glmgam.Rd -552561fc0f9300d1564fb3acb0ab0d84 *man/glmscoretest.Rd +b519af127d6e82fc906146470431867a *man/glmgam.Rd +87ff6241817138438572cd8dd9cc8ff1 *man/glmnbfit.Rd +0c5e6db0bffc33363aad2425cae33e1e *man/glmscoretest.Rd a3395b99cb1d5e2b0a4e65d367e760c9 *man/growthcurve.Rd 0f8f3ac69241aedf8664fbdf9d50f085 *man/hommel.test.Rd f7ff3f35ae819d5ebd1639dcd79a0325 *man/invgauss.Rd @@ -44,14 +45,14 @@ 17ed9e7fcab528fd5ec9631dd8aab8d6 *man/permp.Rd 03485b19a3f73e050a942b597d9addc8 *man/plot.limdil.Rd 8e7ed19b1b289c9ce8fa54a482399d40 *man/power.Rd -4d7a3c0d5daf408d50668b33e7b164e7 *man/qresiduals.Rd +774e2b79c8da8b15447383653ba03da8 *man/qresiduals.Rd 185faf67649fe467b802b8e2fa77fc27 *man/remlscor.Rd e952cfe21c361fc09452f9a621fd52ce *man/remlscorgamma.Rd -396eb7592074f5d79ec1d7e6d9559378 *man/sage.test.Rd +6ad4e994c1bc25987dce8be4f3363299 *man/sage.test.Rd 542b2dac5e344bf19914046c1fb85d57 *man/statmod.Rd 3d0552637454f6b32b1ad801ce684dfe *man/tweedie.Rd d8a58e7e91ba0dbf41668a06d0bc34dc *man/welding.Rd 0848d2cac41cdeea7c861b02f73fafea *src/gaussq2.f f3404859db4cd3831a68abcad9917a9e *src/init.c -f2c61b6611a1f2667184b8e24beba294 *tests/statmod-Tests.R -bf1f58e9bb082bd0f1bf0733c495306a *tests/statmod-Tests.Rout.save +7d9f90d14fa53376fb47c4ffdb79dc86 *tests/statmod-Tests.R +36e792a2c06d66c54cb42ea5623a4ce9 *tests/statmod-Tests.Rout.save diff -Nru r-cran-statmod-1.4.30/NAMESPACE r-cran-statmod-1.4.32/NAMESPACE --- r-cran-statmod-1.4.30/NAMESPACE 2016-11-17 05:11:31.000000000 +0000 +++ r-cran-statmod-1.4.32/NAMESPACE 2019-05-29 07:21:13.000000000 +0000 @@ -7,7 +7,7 @@ importFrom("graphics", "abline", "legend", "lines", "plot", "points") importFrom("stats", "Gamma", "binomial", "chisq.test", "dbinom", - "fisher.test", "fitted", "glm", "glm.fit", "lm.fit", + "fisher.test", "fitted", "glm", "glm.fit", ".lm.fit", "lm.fit", "lm.wfit", "make.link", "median", "model.matrix", "model.response", "model.weights", "p.adjust", "pbeta", "pbinom", "pchisq", "pgamma", "pnorm", "ppois", diff -Nru r-cran-statmod-1.4.30/R/glmgam.R r-cran-statmod-1.4.32/R/glmgam.R --- r-cran-statmod-1.4.30/R/glmgam.R 2012-07-20 05:59:37.000000000 +0000 +++ r-cran-statmod-1.4.32/R/glmgam.R 2019-05-29 06:40:56.000000000 +0000 @@ -1,11 +1,11 @@ # GLMGAM.R -glmgam.fit <- function(X,y,coef.start=NULL,tol=1e-6,maxit=50,trace=FALSE) { +glmgam.fit <- function(X,y,coef.start=NULL,tol=1e-6,maxit=50,trace=FALSE) # Fit gamma generalized linear model with identity link -# by Levenberg damped Fisher scoring +# by Fisher scoring with Levenberg-style damping # Gordon Smyth -# 12 Mar 2003. Last revised 3 November 2010. - +# Created 12 Mar 2003. Last revised 3 November 2010. +{ # check input X <- as.matrix(X) n <- nrow(X) diff -Nru r-cran-statmod-1.4.30/R/glmnb.R r-cran-statmod-1.4.32/R/glmnb.R --- r-cran-statmod-1.4.30/R/glmnb.R 2013-11-20 07:40:58.000000000 +0000 +++ r-cran-statmod-1.4.32/R/glmnb.R 2019-05-29 07:19:20.000000000 +0000 @@ -1,8 +1,8 @@ glmnb.fit <- function(X,y,dispersion,weights=NULL,offset=0,coef.start=NULL,start.method="mean",tol=1e-6,maxit=50,trace=FALSE) # Fit negative binomial generalized linear model with log link -# by Levenberg damped Fisher scoring -# Yunshun Chen and Gordon Smyth -# 2 November 2010. Last modified 20 November 2012. +# by Fisher scoring with Levenberg-style damped +# Gordon Smyth and Yunshun Chen +# Created 2 November 2010. Last modified 29 May 2019. { # Check input values for y y <- as.vector(y) @@ -20,23 +20,51 @@ if(!all(is.finite(X))) stop("All X values must be finite and non-missing") p <- ncol(X) if(p > n) stop("More columns than rows in X") - -# Handle y all zero as special case - if(ymax==0) return(list(coefficients=rep(0,p),fitted.values=rep(0,n),deviance=0,iter=0,convergence="converged")) + if(is.null(colnames(X))) colnames(X) <- paste0("x",1:p) # Check input values for dispersion if(any(dispersion<0)) stop("dispersion values must be non-negative") - phi <- rep(dispersion,length.out=n) + phi <- rep_len(dispersion,n) # Check input values for offset - offset <- rep(offset,length=n) + if(!all(is.finite(offset))) stop("All offset values must be finite and non-missing") + offset <- rep_len(offset,n) # Check input values for weights - if(is.null(weights)) weights <- rep.int(1,n) + if(is.null(weights)) weights <- rep_len(1,n) if(any(weights <= 0)) stop("All weights must be positive") +# Handle y all zero as special case + if(ymax==0) { +# Does X include an intercept term? + if(colnames(X)[1]=="(Intercept)") { + beta <- rep_len(0,p) + names(beta) <- colnames(X) + beta[1] <- -Inf + mu <- rep.int(0,n) + names(mu) <- rownames(X) + return(list(coefficients=beta,fitted.values=mu,deviance=0,iter=0L,convergence="converged")) + } +# Does X span the intercept term, at least closely enough to preserve signs? + One <- rep_len(1,n) + fit <- .lm.fit(X,One) + if(max(abs(fit$residuals)) < 1) { + beta <- -1e10 * fit$coefficients + names(beta) <- colnames(X) + mu <- rep_len(0,n) + names(mu) <- rownames(X) + return(list(coefficients=beta,fitted.values=mu,deviance=0,iter=0L,convergence="converged")) + } +# If X is far from spanning the intercept term, then +# initialize the iteration by trying to cancel out the offsets + if(is.null(coef.start)) { + fit <- lm.wfit(x=X,y=offset,w=weights) + coef.start <- -fit$coefficients + } + } + # Starting values - delta <- min(ymax,1/6) + delta <- 1/6 y1 <- pmax(y,delta) if(is.null(coef.start)) { start.method <- match.arg(start.method,c("log(y)","mean")) @@ -49,7 +77,7 @@ rate <- y/N w <- weights*N/(1+phi*N) beta.mean <- log(sum(w*rate)/sum(w)) - beta <- qr.coef(qr(X),rep(beta.mean,length=n)) + beta <- qr.coef(qr(X),rep_len(beta.mean,n)) mu <- drop(exp(X %*% beta + offset)) } } else { diff -Nru r-cran-statmod-1.4.30/tests/statmod-Tests.R r-cran-statmod-1.4.32/tests/statmod-Tests.R --- r-cran-statmod-1.4.30/tests/statmod-Tests.R 2017-02-27 04:28:42.000000000 +0000 +++ r-cran-statmod-1.4.32/tests/statmod-Tests.R 2019-05-29 07:27:29.000000000 +0000 @@ -16,10 +16,14 @@ glmgam.fit(X=cbind(1,c(1,0.5,0.5,0,0)),y=rchisq(5,df=1)) ### glmnb.fit + y <- rnbinom(5,mu=10,size=10) glmnb.fit(X=cbind(1,c(1,0.5,0.5,0,0)),y=y,dispersion=0.1) glmnb.fit(X=cbind(1,c(1,0.5,0.5,0,0)),y=y,dispersion=runif(6)) glmnb.fit(X=cbind(1,c(1,1,0,0,0)),y=c(0,0,6,2,9),dispersion=0.1) +glmnb.fit(X=cbind(1,c(1,1,0,0,0)),y=c(0,0,0,0,0),dispersion=0.1) +X <- matrix(rnorm(10),5,2) +glmnb.fit(X,y=c(0,0,0,0,0),offset=rnorm(5),dispersion=0.05) ### mixedModel2 diff -Nru r-cran-statmod-1.4.30/tests/statmod-Tests.Rout.save r-cran-statmod-1.4.32/tests/statmod-Tests.Rout.save --- r-cran-statmod-1.4.30/tests/statmod-Tests.Rout.save 2017-02-27 04:30:44.000000000 +0000 +++ r-cran-statmod-1.4.32/tests/statmod-Tests.Rout.save 2019-05-29 07:51:01.000000000 +0000 @@ -1,6 +1,6 @@ -R Under development (unstable) (2017-02-25 r72256) -- "Unsuffered Consequences" -Copyright (C) 2017 The R Foundation for Statistical Computing +R version 3.6.0 (2019-04-26) -- "Planting of a Tree" +Copyright (C) 2019 The R Foundation for Statistical Computing Platform: x86_64-w64-mingw32/x64 (64-bit) R is free software and comes with ABSOLUTELY NO WARRANTY. @@ -82,10 +82,12 @@ > > ### glmnb.fit +> > y <- rnbinom(5,mu=10,size=10) > glmnb.fit(X=cbind(1,c(1,0.5,0.5,0,0)),y=y,dispersion=0.1) $coefficients -[1] 2.3042476 -0.2210662 + x1 x2 + 2.3042476 -0.2210662 $fitted.values [1] 8.029975 8.968465 8.968465 10.016639 10.016639 @@ -101,7 +103,8 @@ > glmnb.fit(X=cbind(1,c(1,0.5,0.5,0,0)),y=y,dispersion=runif(6)) $coefficients -[1] 2.2854591 -0.2049791 + x1 x2 + 2.2854591 -0.2049791 $fitted.values [1] 8.008312 8.872615 8.872615 9.830198 9.830198 @@ -117,7 +120,8 @@ > glmnb.fit(X=cbind(1,c(1,1,0,0,0)),y=c(0,0,6,2,9),dispersion=0.1) $coefficients -[1] 1.734601 -17.510821 + x1 x2 + 1.734601 -17.510821 $fitted.values [1] 1.407586e-07 1.407586e-07 5.666667e+00 5.666667e+00 5.666667e+00 @@ -131,6 +135,41 @@ $convergence [1] "converged" +> glmnb.fit(X=cbind(1,c(1,1,0,0,0)),y=c(0,0,0,0,0),dispersion=0.1) +$coefficients + x1 x2 +-1.000000e+10 2.482534e-06 + +$fitted.values +[1] 0 0 0 0 0 + +$deviance +[1] 0 + +$iter +[1] 0 + +$convergence +[1] "converged" + +> X <- matrix(rnorm(10),5,2) +> glmnb.fit(X,y=c(0,0,0,0,0),offset=rnorm(5),dispersion=0.05) +$coefficients + x1 x2 + 9316725672 -10048340530 + +$fitted.values +[1] 0 0 0 0 0 + +$deviance +[1] 0 + +$iter +[1] 0 + +$convergence +[1] "converged" + > > ### mixedModel2 > @@ -141,28 +180,29 @@ > m$reml.residuals <- m$qr <- NULL > m $varcomp - Residual Block - 0.4541410 -0.1830805 + Residual Block + 2.548669 -0.870409 $se.varcomp -[1] 0.3764600 0.1997047 +[1] 2.543947 1.363837 $coefficients (Intercept) x - -0.1622120 -0.1635664 + 0.1585957 0.5996677 $residuals -[1] 0.18728174 -0.23383906 0.04655731 0.38053667 0.93837501 0.56278279 +[1] 0.002963361 -0.874026581 0.871063220 -1.182699945 1.934015605 +[6] 0.501602395 $fitted.values -[1] 0.42428865 0.43509063 0.06831659 -0.08278783 0.09651970 0.07022039 +[1] -0.2467016 0.0413372 0.1081430 0.8190985 0.1437542 1.0164005 $effects (Intercept) x - -1.8057315 -1.0181717 -0.8221429 0.4683771 1.5047294 0.9167959 + 0.06245044 -0.87448215 1.19519405 -0.26224188 1.29543726 0.90806520 $weights -[1] 11.36623 11.36623 11.36623 2.20196 2.20196 2.20196 +[1] 1.2378524 1.2378524 1.2378524 0.3923617 0.3923617 0.3923617 $rank [1] 2 @@ -174,7 +214,7 @@ [1] 4 $se.coefficients -[1] 0.1331802 0.1606472 +[1] 0.3983904 0.6857404 > > ### mixedModel2Fit @@ -227,14 +267,14 @@ > y <- rnorm(6) > fit <- glm(y~1) > residuals(fit) - 1 2 3 4 5 6 - 0.1570945 -1.0198715 0.6872330 -1.1702352 1.7359615 -0.3901824 + 1 2 3 4 5 6 + 0.68815664 0.33141358 0.07456884 0.39104513 -0.87533184 -0.60985235 > qresiduals(fit) 1 2 3 4 5 6 - 0.1425500 -0.9254473 0.6236059 -1.0618896 1.5752385 -0.3540575 + 1.1222606 0.5404764 0.1216085 0.6377248 -1.4275100 -0.9945603 > qresiduals(fit,dispersion=1) - 1 2 3 4 5 6 - 0.1570945 -1.0198715 0.6872330 -1.1702352 1.7359615 -0.3901824 + 1 2 3 4 5 6 + 0.68815664 0.33141358 0.07456884 0.39104513 -0.87533184 -0.60985235 > > if(require("MASS")) { + fit <- glm(Days~Age,family=negative.binomial(2),data=quine) @@ -244,9 +284,9 @@ + } Loading required package: MASS Min. 1st Qu. Median Mean 3rd Qu. Max. --3.2820 -0.8242 -0.2252 -0.1500 0.7333 3.0359 +-2.9227 -0.8494 -0.2115 -0.1294 0.7212 3.0678 Min. 1st Qu. Median Mean 3rd Qu. Max. --2.90943 -0.52601 -0.02938 -0.01204 0.67881 2.45587 +-3.14845 -0.50446 -0.02932 0.00518 0.67937 2.47162 > > ### gauss.quad > @@ -357,7 +397,7 @@ > qinvgauss(log(p),mean=1.3,dispersion=0.6,lower.tail=FALSE,log.p=TRUE) [1] Inf 9.2602074131 0.9446753861 0.1271035164 0.0000000000 > rinvgauss(5,mean=c(1,NA,3,Inf,1e10),disp=c(2,3,NA,Inf,4)) -[1] 1.122583604 NA NA 0.000000000 36.709796139 +[1] 0.64715825862 NA NA 0.00000000000 0.08624417187 > > ### extra tests done only locally > @@ -368,4 +408,4 @@ > > proc.time() user system elapsed - 0.84 0.03 0.85 + 0.14 0.06 0.20