Improve documentation

Project.toml

@@ -28,7 +28,7 @@ CategoricalArrays = "0.8, 0.9, 0.10, 1"
 DataFrames = "0.22, 1"
 Distributions = "0.16, 0.17, 0.18, 0.19, 0.20, 0.21, 0.22, 0.23, 0.24, 0.25"
 LogExpFunctions = "0.3"
-RDatasets = "0.5, 0.6, 0.7, 0.8"
+RDatasets = "0.8.1"
 Reexport = "0.1, 0.2, 1.0"
 SpecialFunctions = "0.6, 0.7, 0.8, 0.9, 0.10, 1, 2.0"
 Statistics = "1"

docs/make.jl

@@ -5,7 +5,9 @@ makedocs(; format=Documenter.HTML(),
          modules=[GLM],
          pages=["Home" => "index.md",
                 "examples.md",
-                "api.md"],
+                "r-comparison.md",
+                "api.md",
+                "implementation.md"],
          debug=false,
          doctest=true,
          warnonly=[:missing_docs])

docs/src/api.md

@@ -68,6 +68,8 @@ GLM.ftest
 StatsBase.nobs
 StatsBase.nulldeviance
 StatsBase.predict
+StatsBase.stderror
+StatsBase.vcov
 ```
 
 ## Links and methods applied to them

docs/src/implementation.md

@@ -0,0 +1,38 @@
+# Implementation
+
+## Separation of response object and predictor object
+
+The general approach in this code is to separate functionality related
+to the response from that related to the linear predictor. This
+allows for greater generality by mixing and matching different
+subtypes of the abstract type `LinPred` and the abstract type `ModResp`.
+
+A `LinPred` type incorporates the parameter vector and the model
+matrix. The parameter vector is a dense numeric vector but the model
+matrix can be dense or sparse. A `LinPred` type must incorporate
+some form of a decomposition of the weighted model matrix that allows
+for the solution of a system `X'W * X * delta=X'wres` where `W` is a
+diagonal matrix of "X weights", provided as a vector of the square
+roots of the diagonal elements, and `wres` is a weighted residual vector.
+
+Currently there are two dense predictor types, `DensePredQR` and
+`DensePredChol`, and the usual caveats apply. The Cholesky
+version is faster but somewhat less accurate than that QR version.
+The skeleton of a distributed predictor type is in the code
+but not yet fully fleshed out. Because Julia by default uses
+OpenBLAS, which is already multi-threaded on multicore machines, there
+may not be much advantage in using distributed predictor types.
+
+A `ModResp` type must provide methods for the `wtres` and
+`sqrtxwts` generics. Their values are the arguments to the
+`updatebeta` methods of the `LinPred` types. The
+`Float64` value returned by `updatedelta` is the value of the
+convergence criterion.
+
+Similarly, `LinPred` types must provide a method for the
+`linpred` generic. In general `linpred` takes an instance of
+a `LinPred` type and a step factor. Methods that take only an instance
+of a `LinPred` type use a default step factor of 1. The value of
+`linpred` is the argument to the `updatemu` method for
+`ModResp` types. The `updatemu` method returns the updated
+deviance.

docs/src/index.md

@@ -10,12 +10,14 @@ Pkg.add("GLM")
 
 will install this package and its dependencies, which includes the [Distributions package](https://github.com/JuliaStats/Distributions.jl).
 
-The [RDatasets package](https://github.com/johnmyleswhite/RDatasets.jl) is useful for fitting models on standard R datasets to compare the results with those from R.
+The [RDatasets package](https://github.com/JuliaStats/RDatasets.jl) is useful for fitting models on standard R datasets to compare the results with those from R.
 
-## Fitting GLM models
+## Fitting models
+
+Two methods taking different kinds of arguments can be used to fit a model:
+- for linear models: `lm(formula, data)` and `lm(X, y)`;
+- for generalized linear models (GLM): `glm(formula, data, family, link)` and `glm(X, y, family, link)`.
 
-Two methods can be used to fit a Generalized Linear Model (GLM):
-`glm(formula, data, family, link)` and `glm(X, y, family, link)`.
 Their arguments must be:
 
 - `formula`: a [StatsModels.jl `Formula` object](https://juliastats.org/StatsModels.jl/stable/formula/)
@@ -147,6 +149,15 @@ one of the three weights types, and then passed to the `weights` keyword argumen
 Short-hand functions `aweights`, `fweights`, and `pweights` can be used to construct
 `AnalyticWeights`, `FrequencyWeights`, and `ProbabilityWeights`, respectively.
 
+Using analytic weights corresponds to weighted least squares.
+This gives the same results as R and Stata.
+
+Probability weights give the same point estimates as analytic weights, but standard errors
+and p-values are based on a sandwich (heteroskedasticity-robust) estimator.
+This gives the same results as the R `survey` package with a simple survey design
+without strata nor clustering, but differs from Stata with the `pweights` option, which
+adopts the same approach but with a different assumption regarding degrees of freedom.
+
 We illustrate the API with randomly generated data.
 
 ```jldoctest weights
@@ -208,13 +219,10 @@ x            -0.0478667   0.0265532  -1.80    0.0745  -0.100561  0.00482739
 
 ```
 
-!!! warning
-
-In the old API, weights were passed as `AbstractVectors` and were silently treated in
-the internal computation of standard errors and related quantities as `FrequencyWeights`.
-Passing weights as `AbstractVector` is still allowed for backward compatibility, but it
-is deprecated. When weights are passed following the old API, they are now coerced to
-`FrequencyWeights` and a deprecation warning is issued.
+!!! note
+    In the old API before GLM.jl 2.0, weights were passed as `AbstractVectors`
+    and were silently treated in the internal computation of standard errors
+    and related quantities as `FrequencyWeights`.
 
 The type of the weights will affect the variance of the estimated coefficients and the
 quantities involving this variance. The coefficient point estimates will be the same
@@ -282,9 +290,6 @@ Many of the methods provided by this package have names similar to those in [R](
 - `stderror`: standard errors of the coefficients
 - `vcov`: variance-covariance matrix of the coefficient estimates
 
-Note that the canonical link for negative binomial regression is `NegativeBinomialLink`,
-but in practice one typically uses `LogLink`.
-
 ```jldoctest methods
 julia> using GLM, DataFrames, StatsBase
 
@@ -328,43 +333,6 @@ julia> round.(cooksdistance(mdl); digits=8)
  2.5
 ```
 
-## Separation of response object and predictor object
-
-The general approach in this code is to separate functionality related
-to the response from that related to the linear predictor. This
-allows for greater generality by mixing and matching different
-subtypes of the abstract type `LinPred` and the abstract type `ModResp`.
-
-A `LinPred` type incorporates the parameter vector and the model
-matrix. The parameter vector is a dense numeric vector but the model
-matrix can be dense or sparse. A `LinPred` type must incorporate
-some form of a decomposition of the weighted model matrix that allows
-for the solution of a system `X'W * X * delta=X'wres` where `W` is a
-diagonal matrix of "X weights", provided as a vector of the square
-roots of the diagonal elements, and `wres` is a weighted residual vector.
-
-Currently there are two dense predictor types, `DensePredQR` and
-`DensePredChol`, and the usual caveats apply. The Cholesky
-version is faster but somewhat less accurate than that QR version.
-The skeleton of a distributed predictor type is in the code
-but not yet fully fleshed out. Because Julia by default uses
-OpenBLAS, which is already multi-threaded on multicore machines, there
-may not be much advantage in using distributed predictor types.
-
-A `ModResp` type must provide methods for the `wtres` and
-`sqrtxwts` generics. Their values are the arguments to the
-`updatebeta` methods of the `LinPred` types. The
-`Float64` value returned by `updatedelta` is the value of the
-convergence criterion.
-
-Similarly, `LinPred` types must provide a method for the
-`linpred` generic. In general `linpred` takes an instance of
-a `LinPred` type and a step factor. Methods that take only an instance
-of a `LinPred` type use a default step factor of 1. The value of
-`linpred` is the argument to the `updatemu` method for
-`ModResp` types. The `updatemu` method returns the updated
-deviance.
-
 ## Debugging failed fits
 
 In the rare cases when a fit of a generalized linear model fails, it can be useful