From d60a184b205584112496a5ac26b8964623a733f5 Mon Sep 17 00:00:00 2001 From: Mauricio 'Pacha' Vargas Sepulveda Date: Wed, 28 Aug 2024 11:23:58 -0400 Subject: [PATCH] partial progress with srrstatsTODO --- DESCRIPTION | 3 +-- NEWS.md | 1 + R/fixed_effects.R | 6 ++++- R/srr-stats-standards.R | 20 ++++++++-------- tests/testthat/test-felm.R | 13 +++++++++++ tests/testthat/test-fixed-effects.R | 36 +++++++++++++++++++++++++++++ 6 files changed, 66 insertions(+), 13 deletions(-) create mode 100644 tests/testthat/test-fixed-effects.R diff --git a/DESCRIPTION b/DESCRIPTION index cf4c0f8..92abc3f 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -37,10 +37,9 @@ LazyData: true RoxygenNote: 7.3.1 Encoding: UTF-8 NeedsCompilation: yes -LinkingTo: cpp11, armadillo +LinkingTo: cpp11, cpp11armadillo VignetteBuilder: knitr Config/testthat/edition: 3 Remotes: - pachadotdev/armadillo, ropenscilabs/srr Roxygen: list(markdown = TRUE, roclets = c("namespace", "rd", "srr::srr_stats_roclet")) diff --git a/NEWS.md b/NEWS.md index 1c31911..a26a494 100644 --- a/NEWS.md +++ b/NEWS.md @@ -4,6 +4,7 @@ results to R. Previously, there were multiple data copies between R and C++ that added overhead to the computations. * The previous versions returned MX by default, now it has to be specified. +* Adds code to extract the fixed effects with `felm` objects. # capybara 0.5.2 diff --git a/R/fixed_effects.R b/R/fixed_effects.R index 47e4972..81eff79 100644 --- a/R/fixed_effects.R +++ b/R/fixed_effects.R @@ -61,7 +61,11 @@ fixed_effects <- function(object = NULL, alpha_tol = 1.0e-08) { k_list <- get_index_list_(k_vars, data) # Recover fixed effects by alternating the solutions of normal equations - pie <- eta - X %*% beta + if (inherits(object, "feglm")) { + pie <- eta - X %*% beta + } else { + pie <- fitted.values(object) - X %*% beta + } fe_list <- as.list(get_alpha_(pie, k_list, alpha_tol)) # Assign names to the different fixed effects categories diff --git a/R/srr-stats-standards.R b/R/srr-stats-standards.R index 7c8018f..16e50c8 100644 --- a/R/srr-stats-standards.R +++ b/R/srr-stats-standards.R @@ -53,15 +53,14 @@ #' @srrstats {G5.2} *Appropriate error and warning behaviour of all functions should be explicitly demonstrated through tests. In particular,* #' @srrstats {G5.2a} *Every message produced within R code by `stop()`, `warning()`, `message()`, or equivalent should be unique* #' @srrstats {G5.2b} *Explicit tests should demonstrate conditions which trigger every one of those messages, and should compare the result with expected values.* -#' @srrstatsTODO {G5.3} *For functions which are expected to return objects containing no missing (`NA`) or undefined (`NaN`, `Inf`) values, the absence of any such values in return objects should be explicitly tested.* #' @srrstats {G5.4a} *For new methods, it can be difficult to separate out correctness of the method from the correctness of the implementation, as there may not be reference for comparison. In this case, testing may be implemented against simple, trivial cases or against multiple implementations such as an initial R implementation compared with results from a C/C++ implementation.* #' @srrstats {G5.4b} *For new implementations of existing methods, correctness tests should include tests against previous implementations. Such testing may explicitly call those implementations in testing, preferably from fixed-versions of other software, or use stored outputs from those where that is not possible.* #' @srrstats {G5.4c} *Where applicable, stored values may be drawn from published paper outputs when applicable and where code from original implementations is not available* -#' @srrstatsTODO {G5.6} **Parameter recovery tests** *to test that the implementation produce expected results given data with known properties. For instance, a linear regression algorithm should return expected coefficient values for a simulated data set generated from a linear model.* -#' @srrstatsTODO {G5.6a} *Parameter recovery tests should generally be expected to succeed within a defined tolerance rather than recovering exact values.* -#' @srrstatsTODO {G5.6b} *Parameter recovery tests should be run with multiple random seeds when either data simulation or the algorithm contains a random component. (When long-running, such tests may be part of an extended, rather than regular, test suite; see G5.10-4.12, below).* +#' @srrstats {G5.6} **Parameter recovery tests** *to test that the implementation produce expected results given data with known properties. For instance, a linear regression algorithm should return expected coefficient values for a simulated data set generated from a linear model.* +#' @srrstats {G5.6a} *Parameter recovery tests should generally be expected to succeed within a defined tolerance rather than recovering exact values.* +#' @srrstats {G5.6b} *Parameter recovery tests should be run with multiple random seeds when either data simulation or the algorithm contains a random component. (When long-running, such tests may be part of an extended, rather than regular, test suite; see G5.10-4.12, below).* #' @srrstats {G5.7} **Algorithm performance tests** *to test that implementation performs as expected as properties of data change. For instance, a test may show that parameters approach correct estimates within tolerance as data size increases, or that convergence times decrease for higher convergence thresholds.* -#' @srrstatsTODO {G5.8} **Edge condition tests** *to test that these conditions produce expected behaviour such as clear warnings or errors when confronted with data with extreme properties including but not limited to:* +#' @srrstats {G5.8} **Edge condition tests** *to test that these conditions produce expected behaviour such as clear warnings or errors when confronted with data with extreme properties including but not limited to:* #' @srrstats {G5.8a} *Zero-length data* #' @srrstats {G5.8b} *Data of unsupported types (e.g., character or complex numbers in for functions designed only for numeric data)* #' @srrstats {G5.8c} *Data with all-`NA` fields or columns or all identical fields or columns* @@ -70,7 +69,7 @@ #' @srrstatsTODO {G5.9a} *Adding trivial noise (for example, at the scale of `.Machine$double.eps`) to data does not meaningfully change results* #' @srrstats {G5.10} *Extended tests should included and run under a common framework with other tests but be switched on by flags such as as a `_EXTENDED_TESTS="true"` environment variable.* - The extended tests can be then run automatically by GitHub Actions for example by adding the following to the `env` section of the workflow: #' @srrstats {G5.11} *Where extended tests require large data sets or other assets, these should be provided for downloading and fetched as part of the testing workflow.* -#' @srrstatsTODO {G5.12} *Any conditions necessary to run extended tests such as platform requirements, memory, expected runtime, and artefacts produced that may need manual inspection, should be described in developer documentation such as a `CONTRIBUTING.md` or `tests/README.md` file.* +#' @srrstats {G5.12} *Any conditions necessary to run extended tests such as platform requirements, memory, expected runtime, and artefacts produced that may need manual inspection, should be described in developer documentation such as a `CONTRIBUTING.md` or `tests/README.md` file.* #' @srrstats {RE1.0} *Regression Software should enable models to be specified via a formula interface, unless reasons for not doing so are explicitly documented.* #' @srrstats {RE1.1} *Regression Software should document how formula interfaces are converted to matrix representations of input data.* #' @srrstats {RE1.2} *Regression Software should document expected format (types or classes) for inputting predictor variables, including descriptions of types or classes which are not accepted.* @@ -103,7 +102,6 @@ #' @srrstatsTODO {RE4.13} *Predictor variables, and associated "metadata" where applicable.* #' @srrstatsTODO {RE4.14} *Where possible, values should also be provided for extrapolation or forecast *errors*.* #' @srrstatsTODO {RE4.15} *Sufficient documentation and/or testing should be provided to demonstrate that forecast errors, confidence intervals, or equivalent values increase with forecast horizons.* -#' @srrstatsTODO {RE4.16} *Regression Software which models distinct responses for different categorical groups should include the ability to submit new groups to `predict()` methods.* #' @srrstatsTODO {RE4.17} *Model objects returned by Regression Software should implement or appropriately extend a default `print` method which provides an on-screen summary of model (input) parameters and (output) coefficients.* #' @srrstatsTODO {RE4.18} *Regression Software may also implement `summary` methods for model objects, and in particular should implement distinct `summary` methods for any cases in which calculation of summary statistics is computationally non-trivial (for example, for bootstrapped estimates of confidence intervals).* #' @srrstatsTODO {RE5.0} *Scaling relationships between sizes of input data (numbers of observations, with potential extension to numbers of variables/columns) and speed of algorithm.* @@ -113,10 +111,10 @@ #' @srrstatsTODO {RE6.3} *Where a model object is used to generate a forecast (for example, through a `predict()` method), the default `plot` method should provide clear visual distinction between modelled (interpolated) and forecast (extrapolated) values.* #' @srrstatsTODO {RE7.0} *Tests with noiseless, exact relationships between predictor (independent) data.* #' @srrstatsTODO {RE7.0a} In particular, these tests should confirm ability to reject perfectly noiseless input data. -#' @srrstatsTODO {RE7.1} *Tests with noiseless, exact relationships between predictor (independent) and response (dependent) data.* +#' @srrstats {RE7.1} *Tests with noiseless, exact relationships between predictor (independent) and response (dependent) data.* #' @srrstatsTODO {RE7.1a} *In particular, these tests should confirm that model fitting is at least as fast or (preferably) faster than testing with equivalent noisy data (see RE2.4b).* -#' @srrstatsTODO {RE7.2} Demonstrate that output objects retain aspects of input data such as row or case names (see **RE1.3**). -#' @srrstatsTODO {RE7.3} Demonstrate and test expected behaviour when objects returned from regression software are submitted to the accessor methods of **RE4.2**--**RE4.7**. +#' @srrstats {RE7.2} Demonstrate that output objects retain aspects of input data such as row or case names (see **RE1.3**). +#' @srrstats {RE7.3} Demonstrate and test expected behaviour when objects returned from regression software are submitted to the accessor methods of **RE4.2**--**RE4.7**. #' @srrstatsTODO {RE7.4} Extending directly from **RE4.15**, where appropriate, tests should demonstrate and confirm that forecast errors, confidence intervals, or equivalent values increase with forecast horizons. #' @noRd NULL @@ -129,10 +127,12 @@ NULL #' (These comments may also be deleted at any time.) #' @srrstatsNA {G2.14c} *replace missing data with appropriately imputed values* #' @srrstatsNA {G2.16} *All functions should also provide options to handle undefined values (e.g., `NaN`, `Inf` and `-Inf`), including potentially ignoring or removing such values.* +#' @srrstatsNA {G5.3} *For functions which are expected to return objects containing no missing (`NA`) or undefined (`NaN`, `Inf`) values, the absence of any such values in return objects should be explicitly tested.* #' @srrstatsNA {G5.4} **Correctness tests** *to test that statistical algorithms produce expected results to some fixed test data sets (potentially through comparisons using binding frameworks such as [RStata](https://github.com/lbraglia/RStata)).* #' @srrstatsNA {G5.5} *Correctness tests should be run with a fixed random seed* #' @srrstatsNA {G5.9b} *Running under different random seeds or initial conditions does not meaningfully change results* #' @srrstatsNA {G5.11a} *When any downloads of additional data necessary for extended tests fail, the tests themselves should not fail, rather be skipped and implicitly succeed with an appropriate diagnostic message.* #' @srrstatsNA {RE2.2} *Regression Software should provide different options for processing missing values in predictor and response data. For example, it should be possible to fit a model with no missing predictor data in order to generate values for all associated response points, even where submitted response values may be missing.* +#' @srrstatsNA {RE4.16} *Regression Software which models distinct responses for different categorical groups should include the ability to submit new groups to `predict()` methods.* #' @noRd NULL diff --git a/tests/testthat/test-felm.R b/tests/testthat/test-felm.R index d354189..fe4d1e7 100644 --- a/tests/testthat/test-felm.R +++ b/tests/testthat/test-felm.R @@ -50,4 +50,17 @@ test_that("felm works", { m1 <- felm(mpg ~ wt + qsec | cyl + am | carb, mtcars) expect_equal(round(coef(m1), 2), round(coef(m2)[c(2, 3)], 2)) + + set.seed(200100) + d <- data.frame( + y = rnorm(100), + f = factor(sample(1:2, 1000, replace = TRUE)) + ) + d$x <- 2 * y + + fit <- felm(y ~ x | f, data = d) + s1 <- summary(fit) + expect_equal(s1$r.squared, 1) + expect_equal(s1$adj.r.squared, 1) + expect_equal(s1$cm[4], 0) }) diff --git a/tests/testthat/test-fixed-effects.R b/tests/testthat/test-fixed-effects.R new file mode 100644 index 0000000..bedb40a --- /dev/null +++ b/tests/testthat/test-fixed-effects.R @@ -0,0 +1,36 @@ +test_that("fixed_effects is similar to glm", { + # this is also tested in feglm with real data + # this is about SRR {G5.6b}: + # Parameter recovery tests should be run with multiple random seeds when + # either data simulation or the algorithm contains a random component. + + set.seed(200100) + d <- data.frame( + y = rnorm(100), + x = rnorm(100), + f = factor(sample(1:10, 1000, replace = TRUE)) + ) + + fit1 <- glm(y ~ x + f + 0, data = d) + fit2 <- feglm(y ~ x | f, data = d, family = gaussian()) + + c1 <- unname(coef(fit1)[grep("f", names(coef(fit1)))]) + c2 <- unname(drop(fixed_effects(fit2)$f)) + + expect_equal(round(c1 - c2, 3), rep(0, 10)) + + set.seed(100200) + d <- data.frame( + y = rnorm(100), + x = rnorm(100), + f = factor(sample(1:10, 1000, replace = TRUE)) + ) + + fit1 <- lm(y ~ x + f + 0, data = d) + fit2 <- felm(y ~ x | f, data = d) + + c1 <- unname(coef(fit1)[grep("f", names(coef(fit1)))]) + c2 <- unname(drop(fixed_effects(fit2)$f)) + + expect_equal(round(c1 - c2, 3), rep(0, 10)) +})