R Under development (unstable) (2025-06-30 r88369 ucrt) -- "Unsuffered Consequences" Copyright (C) 2025 The R Foundation for Statistical Computing Platform: x86_64-w64-mingw32/x64 R is free software and comes with ABSOLUTELY NO WARRANTY. You are welcome to redistribute it under certain conditions. Type 'license()' or 'licence()' for distribution details. R is a collaborative project with many contributors. Type 'contributors()' for more information and 'citation()' on how to cite R or R packages in publications. Type 'demo()' for some demos, 'help()' for on-line help, or 'help.start()' for an HTML browser interface to help. Type 'q()' to quit R. > library(distfreereg) > > all.equal.distfreereg <- distfreereg:::all.equal.distfreereg > test_dfr_functions <- distfreereg:::test_dfr_functions > > n <- 1e2 > func <- function(X, theta) theta[1] + theta[2]*X[,1] + theta[3]*X[,2] > set.seed(20250516) > Sig <- rWishart(1, df = n, Sigma = diag(n))[,,1] > X <- matrix(rexp(2*n, rate = 1), nrow = n) > theta <- c(2,5,1) > Y <- distfreereg:::f2ftheta(f = func, X)(theta) + + as.vector(distfreereg:::rmvnorm(n = n, reps = 1, mean = rep(0,n), SqrtSigma = distfreereg:::matsqrt(Sig))) > > df_lm <- as.data.frame(cbind(Y, X, rep(1:10, 10))) > names(df_lm) <- c("z", "x", "y", "g") > set.seed(20250516) > wt <- rexp(n) + 1 > set.seed(20250516) > dfr_form_lm <- distfreereg(test_mean = z ~ x + y, data = df_lm, + method_args = list(weights = wt), + method = "lm", verbose = FALSE, + control = list(return_on_error = FALSE)) > > dfr_form_lm Number of observations: 100 Length of EPSP: 100 Monte Carlo simulations: 10000 Estimated parameter values: (Intercept) x y 5.108e+00 6.254e+00 -1.193e+00 Observed statistics: Stat Value Pr(>Value) MCSE KS 4.694e-01 9.325e-01 2.509e-03 CvM 5.971e-02 8.009e-01 3.993e-03 --- `MCSE' is the Monte Carlo standard error of the estimated p-value. > > set.seed(20250516) > dfr_form_lm_no_weights <- distfreereg(test_mean = z ~ x + y, data = df_lm, + method = "lm") Extracting covariance structure from lm object... Calculating Jacobian from 'lm' object... Calculating Jacobian from lm object... Retrieving fitted values from 'lm' object... Calculating the inverse square root of the covariance matrix... Calculating mu... Ordering observations by simplex method... Calculating transformation anchors... Calculating r_tilde... Calculating residuals... All covariate observations are unique; no grouping done... Calculating empirical partial sum process... Calculating observed statistic(s)... Running Monte Carlo simulation... ...1000 of 10000 ...2000 of 10000 ...3000 of 10000 ...4000 of 10000 ...5000 of 10000 ...6000 of 10000 ...7000 of 10000 ...8000 of 10000 ...9000 of 10000 ...10000 of 10000 > > newdata_lm <- data.frame(a = rnorm(10), b = rnorm(10)) > test_dfr_functions(dfr_form_lm, newdata = newdata_lm) (Intercept) x y 5.11 6.25 -1.19 2.5 % 97.5 % (Intercept) 0.937 9.280 x 3.670 8.840 y -3.100 0.715 1 2 3 4 5 6 7 8 9 10 11 12 13 9.64 11.50 10.40 23.20 13.50 7.93 9.06 4.93 10.70 9.09 7.03 4.68 4.95 14 15 16 17 18 19 20 21 22 23 24 25 26 7.46 8.67 9.02 7.67 4.07 12.60 9.08 17.10 10.30 4.65 5.67 5.56 11.00 27 28 29 30 31 32 33 34 35 36 37 38 39 13.60 5.33 6.26 8.43 24.80 26.60 3.53 16.80 16.00 7.48 9.35 6.77 8.83 40 41 42 43 44 45 46 47 48 49 50 51 52 17.10 6.91 5.00 5.15 4.51 7.18 8.51 8.41 8.55 3.94 18.10 7.82 19.40 53 54 55 56 57 58 59 60 61 62 63 64 65 7.04 18.00 27.00 5.36 5.97 7.07 3.18 9.41 13.50 11.70 4.89 4.02 9.57 66 67 68 69 70 71 72 73 74 75 76 77 78 19.50 7.17 15.20 15.60 14.10 7.60 8.53 24.90 13.20 6.11 15.50 5.78 5.66 79 80 81 82 83 84 85 86 87 88 89 90 91 17.90 7.43 6.12 8.34 3.47 6.59 5.59 14.30 6.73 6.56 15.60 16.10 5.47 92 93 94 95 96 97 98 99 100 12.60 16.30 5.60 7.07 9.59 12.30 12.80 19.00 8.90 1 2 3 4 5 6 7 8 9 10 11 12 13 9.64 11.50 10.40 23.20 13.50 7.93 9.06 4.93 10.70 9.09 7.03 4.68 4.95 14 15 16 17 18 19 20 21 22 23 24 25 26 7.46 8.67 9.02 7.67 4.07 12.60 9.08 17.10 10.30 4.65 5.67 5.56 11.00 27 28 29 30 31 32 33 34 35 36 37 38 39 13.60 5.33 6.26 8.43 24.80 26.60 3.53 16.80 16.00 7.48 9.35 6.77 8.83 40 41 42 43 44 45 46 47 48 49 50 51 52 17.10 6.91 5.00 5.15 4.51 7.18 8.51 8.41 8.55 3.94 18.10 7.82 19.40 53 54 55 56 57 58 59 60 61 62 63 64 65 7.04 18.00 27.00 5.36 5.97 7.07 3.18 9.41 13.50 11.70 4.89 4.02 9.57 66 67 68 69 70 71 72 73 74 75 76 77 78 19.50 7.17 15.20 15.60 14.10 7.60 8.53 24.90 13.20 6.11 15.50 5.78 5.66 79 80 81 82 83 84 85 86 87 88 89 90 91 17.90 7.43 6.12 8.34 3.47 6.59 5.59 14.30 6.73 6.56 15.60 16.10 5.47 92 93 94 95 96 97 98 99 100 12.60 16.30 5.60 7.07 9.59 12.30 12.80 19.00 8.90 NULL (Intercept) x y (Intercept) 4.42 -1.76000 -1.07000 x -1.76 1.70000 -0.00825 y -1.07 -0.00825 0.92400 > > > m <- lm(z ~ x + y, data = df_lm, weights = wt) > > set.seed(20250516) > dfr_lm <- distfreereg(test_mean = m, verbose = FALSE, + control = list(return_on_error = FALSE)) > set.seed(20250516) > dfr_lm_verbose <- distfreereg(test_mean = m, + control = list(return_on_error = FALSE), + override = list(J = dfr_lm[["J"]], + fitted_values = dfr_lm[["fitted_values"]])) 'x' and/or 'y' not found in lm object; refitting linear model... Extracting covariance structure from lm object... Using supplied Jacobian... Using supplied fitted values... Calculating the inverse square root of the covariance matrix... Calculating mu... Ordering observations by simplex method... Calculating transformation anchors... Calculating r_tilde... Calculating residuals... All covariate observations are unique; no grouping done... Calculating empirical partial sum process... Calculating observed statistic(s)... Running Monte Carlo simulation... ...1000 of 10000 ...2000 of 10000 ...3000 of 10000 ...4000 of 10000 ...5000 of 10000 ...6000 of 10000 ...7000 of 10000 ...8000 of 10000 ...9000 of 10000 ...10000 of 10000 > > dfr_lm Number of observations: 100 Length of EPSP: 100 Monte Carlo simulations: 10000 Estimated parameter values: (Intercept) x y 5.108e+00 6.254e+00 -1.193e+00 Observed statistics: Stat Value Pr(>Value) MCSE KS 4.694e-01 9.325e-01 2.509e-03 CvM 5.971e-02 8.009e-01 3.993e-03 --- `MCSE' is the Monte Carlo standard error of the estimated p-value. > test_dfr_functions(dfr_lm, newdata = newdata_lm) (Intercept) x y 5.11 6.25 -1.19 2.5 % 97.5 % (Intercept) 0.937 9.280 x 3.670 8.840 y -3.100 0.715 1 2 3 4 5 6 7 8 9 10 11 12 13 9.64 11.50 10.40 23.20 13.50 7.93 9.06 4.93 10.70 9.09 7.03 4.68 4.95 14 15 16 17 18 19 20 21 22 23 24 25 26 7.46 8.67 9.02 7.67 4.07 12.60 9.08 17.10 10.30 4.65 5.67 5.56 11.00 27 28 29 30 31 32 33 34 35 36 37 38 39 13.60 5.33 6.26 8.43 24.80 26.60 3.53 16.80 16.00 7.48 9.35 6.77 8.83 40 41 42 43 44 45 46 47 48 49 50 51 52 17.10 6.91 5.00 5.15 4.51 7.18 8.51 8.41 8.55 3.94 18.10 7.82 19.40 53 54 55 56 57 58 59 60 61 62 63 64 65 7.04 18.00 27.00 5.36 5.97 7.07 3.18 9.41 13.50 11.70 4.89 4.02 9.57 66 67 68 69 70 71 72 73 74 75 76 77 78 19.50 7.17 15.20 15.60 14.10 7.60 8.53 24.90 13.20 6.11 15.50 5.78 5.66 79 80 81 82 83 84 85 86 87 88 89 90 91 17.90 7.43 6.12 8.34 3.47 6.59 5.59 14.30 6.73 6.56 15.60 16.10 5.47 92 93 94 95 96 97 98 99 100 12.60 16.30 5.60 7.07 9.59 12.30 12.80 19.00 8.90 1 2 3 4 5 6 7 8 9 10 11 12 13 9.64 11.50 10.40 23.20 13.50 7.93 9.06 4.93 10.70 9.09 7.03 4.68 4.95 14 15 16 17 18 19 20 21 22 23 24 25 26 7.46 8.67 9.02 7.67 4.07 12.60 9.08 17.10 10.30 4.65 5.67 5.56 11.00 27 28 29 30 31 32 33 34 35 36 37 38 39 13.60 5.33 6.26 8.43 24.80 26.60 3.53 16.80 16.00 7.48 9.35 6.77 8.83 40 41 42 43 44 45 46 47 48 49 50 51 52 17.10 6.91 5.00 5.15 4.51 7.18 8.51 8.41 8.55 3.94 18.10 7.82 19.40 53 54 55 56 57 58 59 60 61 62 63 64 65 7.04 18.00 27.00 5.36 5.97 7.07 3.18 9.41 13.50 11.70 4.89 4.02 9.57 66 67 68 69 70 71 72 73 74 75 76 77 78 19.50 7.17 15.20 15.60 14.10 7.60 8.53 24.90 13.20 6.11 15.50 5.78 5.66 79 80 81 82 83 84 85 86 87 88 89 90 91 17.90 7.43 6.12 8.34 3.47 6.59 5.59 14.30 6.73 6.56 15.60 16.10 5.47 92 93 94 95 96 97 98 99 100 12.60 16.30 5.60 7.07 9.59 12.30 12.80 19.00 8.90 NULL (Intercept) x y (Intercept) 4.42 -1.76000 -1.07000 x -1.76 1.70000 -0.00825 y -1.07 -0.00825 0.92400 > > stopifnot(all.equal(dfr_lm, dfr_form_lm)) > stopifnot(all.equal(dfr_lm, dfr_lm_verbose)) > > cdfr_form_lm <- asymptotics(dfr_form_lm, reps = 5) Running simulation... Repetition 1 of 5 Repetition 2 of 5 Repetition 3 of 5 Repetition 4 of 5 Repetition 5 of 5 > cdfr_lm <- asymptotics(dfr_lm, reps = 5) Running simulation... Repetition 1 of 5 Repetition 2 of 5 Repetition 3 of 5 Repetition 4 of 5 Repetition 5 of 5 > > signif(rejection(cdfr_form_lm, alpha = c(0.1, 0.5))[,2:3], digits = 3) alpha rate 1 0.1 0.0 2 0.5 0.2 3 0.1 0.0 4 0.5 0.6 > signif(rejection(cdfr_lm, alpha = c(0.1, 0.5))[,2:3], digits = 3) alpha rate 1 0.1 0.0 2 0.5 0.4 3 0.1 0.0 4 0.5 0.4 > > > dfr_lm_but_not_lm <- dfr_lm > class(dfr_lm_but_not_lm[["test_mean"]]) <- "wrong" > tryCatch(asymptotics(dfr_lm_but_not_lm, reps = 5), + error = function(e) warning(e)) Warning message: In FUN(X[[i]], ...) : Invalid input to get_n_single(): wrong > > vcov(dfr_lm, jacobian_args = list("ignored")) (Intercept) x y (Intercept) 4.415110 -1.763615201 -1.074803213 x -1.763615 1.697822450 -0.008247111 y -1.074803 -0.008247111 0.924421379 Warning message: In vcov.distfreereg(dfr_lm, jacobian_args = list("ignored")) : jacobian_args and/or hessian_args ignored when 'test_mean' is not a function > > # Orderings > > set.seed(20250516) > dfr_lm_asis <- update(dfr_lm, ordering = "asis") > set.seed(20250516) > dfr_form_lm_asis <- update(dfr_form_lm, ordering = "asis") > stopifnot(all.equal(dfr_lm_asis, dfr_form_lm_asis)) > > set.seed(20250516) > dfr_lm_optimal <- update(dfr_lm, ordering = "optimal") > set.seed(20250516) > dfr_form_lm_optimal <- update(dfr_form_lm, ordering = "optimal") > stopifnot(all.equal(dfr_lm_optimal, dfr_form_lm_optimal)) > > set.seed(20250516) > dfr_lm_natural <- update(dfr_lm, ordering = "natural") > set.seed(20250516) > dfr_form_lm_natural <- update(dfr_form_lm, ordering = "natural") > stopifnot(all.equal(dfr_lm_natural, dfr_form_lm_natural)) > > dfr_lm_natural_grouped <- update(dfr_lm, ordering = "natural", group = TRUE, + verbose = TRUE) 'x' and/or 'y' not found in lm object; refitting linear model... Extracting covariance structure from lm object... Calculating Jacobian from 'lm' object... Calculating Jacobian from lm object... Retrieving fitted values from 'lm' object... Calculating the inverse square root of the covariance matrix... Calculating mu... Ordering observations by natural order... Calculating transformation anchors... Calculating r_tilde... Calculating residuals... All covariate observations are unique; no grouping done... Calculating empirical partial sum process... Calculating observed statistic(s)... Running Monte Carlo simulation... ...1000 of 10000 ...2000 of 10000 ...3000 of 10000 ...4000 of 10000 ...5000 of 10000 ...6000 of 10000 ...7000 of 10000 ...8000 of 10000 ...9000 of 10000 ...10000 of 10000 > > set.seed(20250516) > dfr_lm_g_character <- update(dfr_lm, ordering = list("g")) > > dfr_lm_g_character Number of observations: 100 Length of EPSP: 10 Monte Carlo simulations: 10000 Estimated parameter values: (Intercept) x y 5.108e+00 6.254e+00 -1.193e+00 Observed statistics: Stat Value Pr(>Value) MCSE KS 4.378e-01 2.884e-01 4.530e-03 CvM 5.995e-02 1.849e-01 3.882e-03 --- `MCSE' is the Monte Carlo standard error of the estimated p-value. > > set.seed(20250516) > dfr_form_lm_g <- update(dfr_form_lm, ordering = list("g")) > stopifnot(all.equal(dfr_lm_g_character, dfr_form_lm_g)) > > df_lm[dfr_lm_g_character[["res_order"]],][["g"]] [1] 1 2 3 4 5 6 7 8 9 10 1 2 3 4 5 6 7 8 9 10 1 2 3 4 5 [26] 6 7 8 9 10 1 2 3 4 5 6 7 8 9 10 1 2 3 4 5 6 7 8 9 10 [51] 1 2 3 4 5 6 7 8 9 10 1 2 3 4 5 6 7 8 9 10 1 2 3 4 5 [76] 6 7 8 9 10 1 2 3 4 5 6 7 8 9 10 1 2 3 4 5 6 7 8 9 10 > df_lm[dfr_form_lm_g[["res_order"]],][["g"]] [1] 1 2 3 4 5 6 7 8 9 10 1 2 3 4 5 6 7 8 9 10 1 2 3 4 5 [26] 6 7 8 9 10 1 2 3 4 5 6 7 8 9 10 1 2 3 4 5 6 7 8 9 10 [51] 1 2 3 4 5 6 7 8 9 10 1 2 3 4 5 6 7 8 9 10 1 2 3 4 5 [76] 6 7 8 9 10 1 2 3 4 5 6 7 8 9 10 1 2 3 4 5 6 7 8 9 10 > > set.seed(20250516) > dfr_lm_g_character_grouped <- update(dfr_lm_g_character, group = TRUE) > set.seed(20250516) > dfr_form_lm_g_grouped <- update(dfr_form_lm_g, group = TRUE) > stopifnot(all.equal(dfr_lm_g_character_grouped, dfr_form_lm_g_grouped)) > > > > ### Failures > > tryCatch(update(dfr_lm, theta_init = 1), error = function(e) warning(e)) Warning message: In validate_extra_arg_list(extra_arg_list, "distfreereg.lm()") : Unused arguments passed to distfreereg.lm() with name(s) theta_init > tryCatch(update(dfr_lm, ordering = list("h")), error = function(e) warning(e)) Warning message: In validate_order_columns(control[["data"]], ordering) : Ordering columns not found in data: h > tryCatch(update(dfr_lm, ordering = list(1, "g")), error = function(e) warning(e)) Warning message: In validate_order_columns(control[["data"]], ordering) : Ordering columns not found in data: 1 > tryCatch(update(dfr_lm, ordering = c(1)), error = function(e) warning(e)) Warning message: In strict_match(ordering, c("asis", "optimal", "simplex", "natural")) : ordering must be exactly one of 'asis', 'optimal', 'simplex', 'natural'; it cannot be '1' > tryCatch(update(dfr_lm, ordering = list(1:10)), error = function(e) warning(e)) Warning message: In validate_order_columns(control[["data"]], ordering) : Ordering column specification must be a list of character elements. > > dfr_lm_fail <- update(dfr_lm, override = list(r = matrix(rnorm(length(dfr_lm[["J"]])), nrow = n)), + control = list(return_on_error = TRUE)) Warning message: Error encountered, partial results returned: Error in doTryCatch(return(expr), name, parentenv, handler): Columns of r are not orthogonal > > tryCatch(update(dfr_form_lm, test_mean = z + y ~ x), error = function(e) warning(e)) Warning message: In validate_single_response_term(form = test_mean) : model formula must have a single response term > > tryCatch(update(dfr_form_lm, group = "hello"), error = function(e) warning(e)) Warning message: In validate_args_distfreereg_default(Y = Y, X = X, covariance = covariance, : 'group' must be TRUE or FALSE > > > proc.time() user system elapsed 13.37 0.79 14.17