This function estimates maximum likelihood models with any number of fixed-effects.

femlm(
fml,
data,
family = c("poisson", "negbin", "logit", "gaussian"),
start = 0,
fixef,
fixef.rm = "perfect",
offset,
subset,
split,
fsplit,
cluster,
se,
dof,
panel.id,
fixef.tol = 1e-05,
fixef.iter = 10000,
lean = FALSE,
verbose = 0,
warn = TRUE,
notes = getFixest_notes(),
theta.init,
combine.quick,
mem.clean = FALSE,
only.env = FALSE,
env,
...
)

fenegbin(
fml,
data,
theta.init,
start = 0,
fixef,
fixef.rm = "perfect",
offset,
subset,
split,
fsplit,
cluster,
se,
dof,
panel.id,
fixef.tol = 1e-05,
fixef.iter = 10000,
lean = FALSE,
verbose = 0,
warn = TRUE,
notes = getFixest_notes(),
combine.quick,
mem.clean = FALSE,
only.env = FALSE,
env,
...
)

## Arguments

fml A formula representing the relation to be estimated. For example: fml = z~x+y. To include fixed-effects, insert them in this formula using a pipe: e.g. fml = z~x+y|fixef_1+fixef_2. Multiple estimations can be performed at once: for multiple dep. vars, wrap them in c(): ex c(y1, y2). For multiple indep. vars, use the stepwise functions: ex x1 + csw(x2, x3). The formula fml = c(y1, y2) ~ x1 + cw0(x2, x3) leads to 6 estimation, see details. A data.frame containing the necessary variables to run the model. The variables of the non-linear right hand side of the formula are identified with this data.frame names. Can also be a matrix. Character scalar. It should provide the family. The possible values are "poisson" (Poisson model with log-link, the default), "negbin" (Negative Binomial model with log-link), "logit" (LOGIT model with log-link), "gaussian" (Gaussian model). Starting values for the coefficients. Can be: i) a numeric of length 1 (e.g. start = 0, the default), ii) a numeric vector of the exact same length as the number of variables, or iii) a named vector of any length (the names will be used to initialize the appropriate coefficients). Character vector. The names of variables to be used as fixed-effects. These variables should contain the identifier of each observation (e.g., think of it as a panel identifier). Note that the recommended way to include fixed-effects is to insert them directly in the formula. Can be equal to "perfect" (default), "singleton", "both" or "none". Controls which observations are to be removed. If "perfect", then observations having a fixed-effect with perfect fit (e.g. only 0 outcomes in Poisson estimations) will be removed. If "singleton", all observations for which a fixed-effect appears only once will be removed. The meaning of "both" and "none" is direct. A formula or a numeric vector. An offset can be added to the estimation. If equal to a formula, it should be of the form (for example) ~0.5*x**2. This offset is linearly added to the elements of the main formula 'fml'. A vector (logical or numeric) or a one-sided formula. If provided, then the estimation will be performed only on the observations defined by this argument. A one sided formula representing a variable (eg split = ~var) or a vector. If provided, the sample is split according to the variable and one estimation is performed for each value of that variable. If you also want to include the estimation for the full sample, use the argument fsplit instead. A one sided formula representing a variable (eg split = ~var) or a vector. If provided, the sample is split according to the variable and one estimation is performed for each value of that variable. This argument is the same as split but also includes the full sample as the first estimation. Tells how to cluster the standard-errors (if clustering is requested). Can be either a list of vectors, a character vector of variable names, a formula or an integer vector. Assume we want to perform 2-way clustering over var1 and var2 contained in the data.frame base used for the estimation. All the following cluster arguments are valid and do the same thing: cluster = base[, c("var1", "var2")], cluster = c("var1", "var2"), cluster = ~var1+var2. If the two variables were used as clusters in the estimation, you could further use cluster = 1:2 or leave it blank with se = "twoway" (assuming var1 [resp. var2] was the 1st [res. 2nd] cluster). You can interact two variables using ^ with the following syntax: cluster = ~var1^var2 or cluster = "var1^var2". Character scalar. Which kind of standard error should be computed: “standard”, “hetero”, “cluster”, “twoway”, “threeway” or “fourway”? By default if there are clusters in the estimation: se = "cluster", otherwise se = "standard". Note that this argument can be implicitly deduced from the argument cluster. An object of class dof.type obtained with the function dof. Represents how the degree of freedom correction should be done.You must use the function dof for this argument. The arguments and defaults of the function dof are: adj = TRUE, fixef.K="nested", cluster.adj = TRUE, cluster.df = "conventional", t.df = "conventional", fixef.force_exact=FALSE). See the help of the function dof for details. The panel identifiers. Can either be: i) a one sided formula (e.g. panel.id = ~id+time), ii) a character vector of length 2 (e.g. panel.id=c('id', 'time'), or iii) a character scalar of two variables separated by a comma (e.g. panel.id='id,time'). Note that you can combine variables with ^ only inside formulas (see the dedicated section in feols). Precision used to obtain the fixed-effects. Defaults to 1e-5. It corresponds to the maximum absolute difference allowed between two coefficients of successive iterations. Argument fixef.tol cannot be lower than 10000*.Machine\$double.eps. Note that this parameter is dynamically controlled by the algorithm. Maximum number of iterations in fixed-effects algorithm (only in use for 2+ fixed-effects). Default is 10000. The number of threads. Can be: a) an integer lower than, or equal to, the maximum number of threads; b) 0: meaning all available threads will be used; c) a number strictly between 0 and 1 which represents the fraction of all threads to use. The default is to use 50% of all threads. You can set permanently the number of threads used within this package using the function setFixest_nthreads. Logical, default is FALSE. If TRUE then all large objects are removed from the returned result: this will save memory but will block the possibility to use many methods. It is recommended to use the arguments se or cluster to obtain the appropriate standard-errors at estimation time, since obtaining different SEs won't be possible afterwards. Integer, default is 0. It represents the level of information that should be reported during the optimisation process. If verbose=0: nothing is reported. If verbose=1: the value of the coefficients and the likelihood are reported. If verbose=2: 1 + information on the computing time of the null model, the fixed-effects coefficients and the hessian are reported. Logical, default is TRUE. Whether warnings should be displayed (concerns warnings relating to convergence state). Logical. By default, two notes are displayed: when NAs are removed (to show additional information) and when some observations are removed because of only 0 (or 0/1) outcomes in a fixed-effect setup (in Poisson/Neg. Bin./Logit models). To avoid displaying these messages, you can set notes = FALSE. You can remove these messages permanently by using setFixest_notes(FALSE). Positive numeric scalar. The starting value of the dispersion parameter if family="negbin". By default, the algorithm uses as a starting value the theta obtained from the model with only the intercept. Logical. When you combine different variables to transform them into a single fixed-effects you can do e.g. y ~ x | paste(var1, var2). The algorithm provides a shorthand to do the same operation: y ~ x | var1^var2. Because pasting variables is a costly operation, the internal algorithm may use a numerical trick to hasten the process. The cost of doing so is that you lose the labels. If you are interested in getting the value of the fixed-effects coefficients after the estimation, you should use combine.quick = FALSE. By default it is equal to FALSE if the number of observations is lower than 50,000, and to TRUE otherwise. Logical, default is FALSE. Only to be used if the data set is large compared to the available RAM. If TRUE then intermediary objects are removed as much as possible and gc is run before each substantial C++ section in the internal code to avoid memory issues. (Advanced users.) Logical, default is FALSE. If TRUE, then only the environment used to make the estimation is returned. (Advanced users.) A fixest environment created by a fixest estimation with only.env = TRUE. Default is missing. If provided, the data from this environment will be used to perform the estimation. Not currently used.

## Value

A fixest object. Note that fixest objects contain many elements and most of them are for internal use, they are presented here only for information. To access them, it is safer to use the user-level methods (e.g. vcov.fixest, resid.fixest, etc) or functions (like for instance fitstat to access any fit statistic).

nobs

The number of observations.

fml

The linear formula of the call.

call

The call of the function.

method

The method used to estimate the model.

family

The family used to estimate the model.

fml_all

A list containing different parts of the formula. Always contain the linear formula. Then, if relevant: fixef: the fixed-effects; NL: the non linear part of the formula.

nparams

The number of parameters of the model.

fixef_vars

The names of each fixed-effect dimension.

fixef_id

The list (of length the number of fixed-effects) of the fixed-effects identifiers for each observation.

fixef_sizes

The size of each fixed-effect (i.e. the number of unique identifierfor each fixed-effect dimension).

convStatus

Logical, convergence status.

message

The convergence message from the optimization procedures.

obsRemoved

(When relevant.) In the case there were fixed-effects and some observations were removed because of only 0/1 outcome within a fixed-effect, it gives the row numbers of the observations that were removed. Also reports the NA observations that were removed.

fixef_removed

(When relevant.) In the case there were fixed-effects and some observations were removed because of only 0/1 outcome within a fixed-effect, it gives the list (for each fixed-effect dimension) of the fixed-effect identifiers that were removed.

coefficients

The named vector of estimated coefficients.

coeftable

The table of the coefficients with their standard errors, z-values and p-values.

loglik

The log-likelihood.

iterations

Number of iterations of the algorithm.

ll_null

Log-likelihood of the null model (i.e. with the intercept only).

ll_fe_only

Log-likelihood of the model with only the fixed-effects.

ssr_null

Sum of the squared residuals of the null model (containing only with the intercept).

pseudo_r2

fitted.values

The fitted values are the expected value of the dependent variable for the fitted model: that is $$E(Y|X)$$.

residuals

The residuals (y minus the fitted values).

sq.cor

Squared correlation between the dependent variable and the expected predictor (i.e. fitted.values) obtained by the estimation.

hessian

The Hessian of the parameters.

cov.unscaled

The variance-covariance matrix of the parameters.

se

The standard-error of the parameters.

scores

The matrix of the scores (first derivative for each observation).

residuals

The difference between the dependent variable and the expected predictor.

sumFE

The sum of the fixed-effects coefficients for each observation.

offset

(When relevant.) The offset formula.

weights

(When relevant.) The weights formula.

## Details

Note that the functions feglm and femlm provide the same results when using the same families but differ in that the latter is a direct maximum likelihood optimization (so the two can really have different convergence rates).

## Combining the fixed-effects

You can combine two variables to make it a new fixed-effect using ^. The syntax is as follows: fe_1^fe_2. Here you created a new variable which is the combination of the two variables fe_1 and fe_2. This is identical to doing paste0(fe_1, "_", fe_2) but more convenient.

Note that pasting is a costly operation, especially for large data sets. Thus, the internal algorithm uses a numerical trick which is fast, but the drawback is that the identity of each observation is lost (i.e. they are now equal to a meaningless number instead of being equal to paste0(fe_1, "_", fe_2)). These “identities” are useful only if you're interested in the value of the fixed-effects (that you can extract with fixef.fixest). If you're only interested in coefficients of the variables, it doesn't matter. Anyway, you can use combine.quick = FALSE to tell the internal algorithm to use paste instead of the numerical trick. By default, the numerical trick is performed only for large data sets.

## Lagging variables

To use leads/lags of variables in the estimation, you can: i) either provide the argument panel.id, ii) either set your data set as a panel with the function panel. Doing either of the two will give you acceess to the lagging functions l, f and d.

You can provide several leads/lags/differences at once: e.g. if your formula is equal to f(y) ~ l(x, -1:1), it means that the dependent variable is equal to the lead of y, and you will have as explanatory variables the lead of x1, x1 and the lag of x1. See the examples in function l for more details.

## Interactions

You can interact a numeric variable with a "factor-like" variable by using interact(var, fe, ref), where fe is the variable to be interacted with and the argument ref is a value of fe taken as a reference (optional). Instead of using the function interact, you can use the alias i(var, fe, ref).

Using this specific way to create interactions leads to a different display of the interacted values in etable and offers a special representation of the interacted coefficients in the function coefplot. See examples.

It is important to note that *if you do not care about the standard-errors of the interactions*, then you can add interactions in the fixed-effects part of the formula (using the syntax fe[[var]], as explained in the section “Varying slopes”).

The function interact has in fact more arguments, please see details in its associated help page.

## On standard-errors

Standard-errors can be computed in different ways, you can use the arguments se and dof in summary.fixest to define how to compute them. By default, in the presence of fixed-effects, standard-errors are automatically clustered.

The following vignette: On standard-errors describes in details how the standard-errors are computed in fixest and how you can replicate standard-errors from other software.

You can use the functions setFixest_se and setFixest_dof to permanently set the way the standard-errors are computed.

## Multiple estimations

Multiple estimations can be performed at once, they just have to be specified in the formula. Multiple estimations yield a fixest_multi object which is ‘kind of’ a list of all the results but includes specific methods to access the results in a handy way.

To include mutliple dependent variables, wrap them in c() (list() also works). For instance fml = c(y1, y2) ~ x1 would estimate the model fml = y1 ~ x1 and then the model fml = y2 ~ x1.

To include multiple independent variables, you need to use the stepwise functions. There are 4 stepwise functions associated to 4 short aliases. These are a) stepwise, stepwise0, cstepwise, cstepwise0, and b) sw, sw0, csw, csw0. Let's explain that. Assume you have the following formula: fml = y ~ x1 + sw(x2, x3). The stepwise function sw will estimate the following two models: y ~ x1 + x2 and y ~ x1 + x3. That is, each element in sw() is sequentially, and separately, added to the formula. Would have you used sw0 in lieu of sw, then the model y ~ x1 would also have been estimated. The 0 in the name means that the model wihtout any stepwise element also needs to be estimated. Finally, the prefix c means cumulative: each stepwise element is added to the next. That is, fml = y ~ x1 + csw(x2, x3) would lead to the following models y ~ x1 + x2 and y ~ x1 + x2 + x3. The 0 has the same meaning and would also lead to the model without the stepwise elements to be estimated: in other words, fml = y ~ x1 + csw0(x2, x3) leads to the following three models: y ~ x1, y ~ x1 + x2 and y ~ x1 + x2 + x3.

Multiple independent variables can be combined with multiple dependent variables, as in fml = c(y1, y2) ~ cw(x1, x2, x3) which would lead to 6 estimations. Multiple estimations can also be combined to split samples (with the arguments split, fsplit).

Fixed-effects cannot be included in a stepwise fashion: they are there or not and stay the same for all estimations.

A note on performance. The feature of multiple estimations has been highly optimized for feols, in particular in the presence of fixed-effects. It is faster to estimate multiple models using the formula rather than with a loop. For non-feols models using the formula is roughly similar to using a loop performance-wise.

## References

Berge, Laurent, 2018, "Efficient estimation of maximum likelihood models with multiple fixed-effects: the R package FENmlm." CREA Discussion Papers, 13 (https://wwwen.uni.lu/content/download/110162/1299525/file/2018_13).

For models with multiple fixed-effects:

Gaure, Simen, 2013, "OLS with multiple high dimensional category variables", Computational Statistics & Data Analysis 66 pp. 8--18

On the unconditionnal Negative Binomial model:

Allison, Paul D and Waterman, Richard P, 2002, "Fixed-Effects Negative Binomial Regression Models", Sociological Methodology 32(1) pp. 247--265

See also summary.fixest to see the results with the appropriate standard-errors, fixef.fixest to extract the fixed-effects coefficients, and the function etable to visualize the results of multiple estimations. And other estimation methods: feols, feglm, fepois, feNmlm.

Laurent Berge

## Examples



# We estimate the effect of distance on trade => we account for 3 fixed-effects
# 1) Poisson estimation
est_pois = femlm(Euros ~ log(dist_km) | Origin + Destination + Product, trade)

# 2) Log-Log Gaussian estimation (with same FEs)
est_gaus = update(est_pois, log(Euros+1) ~ ., family = "gaussian")

# Comparison of the results using the function etable
etable(est_pois, est_gaus)
#>                           est_pois           est_gaus
#> Dependent Var.:              Euros       log(Euros+1)
#>
#> log(dist_km)    -1.528*** (0.1156) -2.168*** (0.1542)
#> Fixed-Effects:  ------------------ ------------------
#> Origin                         Yes                Yes
#> Destination                    Yes                Yes
#> Product                        Yes                Yes
#> _______________ __________________ __________________
#> Family                     Poisson           Gaussian
#> S.E.: Clustered         by: Origin         by: Origin
#> Observations                38,325             38,325
#> Squared Cor.               0.60377            0.70434
#> Pseudo R2                  0.76039            0.23563
#> BIC                       1.43e+12          152,008.4# Now using two way clustered standard-errors
etable(est_pois, est_gaus, se = "twoway")
#>                           est_pois           est_gaus
#> Dependent Var.:              Euros       log(Euros+1)
#>
#> log(dist_km)    -1.528*** (0.1307) -2.168*** (0.1712)
#> Fixed-Effects:  ------------------ ------------------
#> Origin                         Yes                Yes
#> Destination                    Yes                Yes
#> Product                        Yes                Yes
#> _______________ __________________ __________________
#> Family                     Poisson           Gaussian
#> S.E.: Clustered  by: Orig. & Dest.  by: Orig. & Dest.
#> Observations                38,325             38,325
#> Squared Cor.               0.60377            0.70434
#> Pseudo R2                  0.76039            0.23563
#> BIC                       1.43e+12          152,008.4
# Comparing different types of standard errors
sum_hetero   = summary(est_pois, se = "hetero")
sum_oneway   = summary(est_pois, se = "cluster")
sum_twoway   = summary(est_pois, se = "twoway")
sum_threeway = summary(est_pois, se = "threeway")

etable(sum_hetero, sum_oneway, sum_twoway, sum_threeway)
#>                         sum_hetero         sum_oneway         sum_twoway
#> Dependent Var.:              Euros              Euros              Euros
#>
#> log(dist_km)    -1.528*** (0.0220) -1.528*** (0.1156) -1.528*** (0.1307)
#> Fixed-Effects:  ------------------ ------------------ ------------------
#> Origin                         Yes                Yes                Yes
#> Destination                    Yes                Yes                Yes
#> Product                        Yes                Yes                Yes
#> _______________ __________________ __________________ __________________
#> S.E. type       Heteroskedas.-rob.         by: Origin  by: Orig. & Dest.
#> Observations                38,325             38,325             38,325
#> Squared Cor.               0.60377            0.60377            0.60377
#> Pseudo R2                  0.76039            0.76039            0.76039
#> BIC                       1.43e+12           1.43e+12           1.43e+12
#>                       sum_threeway
#> Dependent Var.:              Euros
#>
#> log(dist_km)    -1.528*** (0.1340)
#> Fixed-Effects:  ------------------
#> Origin                         Yes
#> Destination                    Yes
#> Product                        Yes
#> _______________ __________________
#> S.E. type       by: Orig. & Dest..
#> Observations                38,325
#> Squared Cor.               0.60377
#> Pseudo R2                  0.76039
#> BIC                       1.43e+12

#
# Multiple estimations:
#

# 6 estimations
est_mult = femlm(c(Ozone, Solar.R) ~ Wind + Temp + csw0(Wind:Temp, Day), airquality)

# We can display the results for the first lhs:
etable(est_mult[lhs = 1])
#>                             model 1            model 2            model 3
#> Dependent Var.:               Ozone              Ozone              Ozone
#>
#> (Intercept)       0.5334** (0.1943)   0.8890* (0.4430)   0.8028. (0.4441)
#> Wind            -0.0761*** (0.0052) -0.1151** (0.0439) -0.1152** (0.0440)
#> Temp             0.0483*** (0.0020) 0.0441*** (0.0051) 0.0444*** (0.0051)
#> Wind x Temp                            0.0005 (0.0005)    0.0005 (0.0005)
#> Day                                                    0.0051*** (0.0015)
#> _______________ ___________________ __________________ __________________
#> S.E. type                  Standard           Standard           Standard
#> Observations                    116                116                116
#> Squared Cor.                0.62615            0.63129            0.63379
#> Pseudo R2                   0.53377            0.53402            0.53759
#> BIC                         1,564.7            1,568.7            1,561.5
# And now the second (access can be made by name)
etable(est_mult[lhs = "Solar.R"])
#>                            model 1            model 2             model 3
#> Dependent Var.:            Solar.R            Solar.R             Solar.R
#>
#> (Intercept)      3.777*** (0.0714)  4.343*** (0.1689)   4.529*** (0.1708)
#> Wind            0.0119*** (0.0020) -0.0415** (0.0146)  -0.0433** (0.0147)
#> Temp            0.0169*** (0.0008) 0.0098*** (0.0021)  0.0088*** (0.0021)
#> Wind x Temp                        0.0007*** (0.0002)  0.0007*** (0.0002)
#> Day                                                   -0.0061*** (0.0007)
#> _______________ __________________ __________________ ___________________
#> S.E. type                 Standard           Standard            Standard
#> Observations                   146                146                 146
#> Squared Cor.               0.08076            0.08166             0.08983
#> Pseudo R2                  0.06045            0.06207             0.07100
#> BIC                        8,192.6            8,183.5             8,110.7
# Now we focus on the two last right hand sides
# (note that .N can be used to specify the last item)
etable(est_mult[rhs = 2:.N])
#>                            model 1            model 2            model 3
#> Dependent Var.:              Ozone            Solar.R              Ozone
#>
#> (Intercept)       0.8890* (0.4430)  4.343*** (0.1689)   0.8028. (0.4441)
#> Wind            -0.1151** (0.0439) -0.0415** (0.0146) -0.1152** (0.0440)
#> Temp            0.0441*** (0.0051) 0.0098*** (0.0021) 0.0444*** (0.0051)
#> Wind x Temp        0.0005 (0.0005) 0.0007*** (0.0002)    0.0005 (0.0005)
#> Day                                                   0.0051*** (0.0015)
#> _______________ __________________ __________________ __________________
#> S.E. type                 Standard           Standard           Standard
#> Observations                   116                146                116
#> Squared Cor.               0.63129            0.08166            0.63379
#> Pseudo R2                  0.53402            0.06207            0.53759
#> BIC                        1,568.7            8,183.5            1,561.5
#>                             model 4
#> Dependent Var.:             Solar.R
#>
#> (Intercept)       4.529*** (0.1708)
#> Wind             -0.0433** (0.0147)
#> Temp             0.0088*** (0.0021)
#> Wind x Temp      0.0007*** (0.0002)
#> Day             -0.0061*** (0.0007)
#> _______________ ___________________
#> S.E. type                  Standard
#> Observations                    146
#> Squared Cor.                0.08983
#> Pseudo R2                   0.07100
#> BIC                         8,110.7
# Combining with split
est_split = fepois(c(Ozone, Solar.R) ~ sw(poly(Wind, 2), poly(Temp, 2)),
airquality, split = ~ Month)

# You can display everything at once with the print method
est_split
#> Standard-errors: Standard
#>
#>
#> # SAMPLE: 5
#>
#>
#> ### Dep. var.: Ozone
#>
#> Expl. vars.: poly(Wind, 2)
#>                Estimate Std. Error  z value  Pr(>|z|)
#> (Intercept)      3.3040   0.042415  77.8970 < 2.2e-16 ***
#> poly(Wind, 2)1  -6.4147   0.617300 -10.3910 < 2.2e-16 ***
#> poly(Wind, 2)2   3.0758   0.570194   5.3944  6.88e-08 ***
#> ---
#> Expl. vars.: poly(Temp, 2)
#>                Estimate Std. Error z value  Pr(>|z|)
#> (Intercept)      4.0582   0.109394 37.0970 < 2.2e-16 ***
#> poly(Temp, 2)1  11.2850   1.448100  7.7926  6.56e-15 ***
#> poly(Temp, 2)2   2.5057   1.171800  2.1383  0.032492 *
#>
#> ### Dep. var.: Solar.R
#>
#> Expl. vars.: poly(Wind, 2)
#>                 Estimate Std. Error  z value  Pr(>|z|)
#> (Intercept)     5.207800   0.018194 286.2400 < 2.2e-16 ***
#> poly(Wind, 2)1 -0.569626   0.271696  -2.0966  0.036033 *
#> poly(Wind, 2)2 -1.834700   0.250895  -7.3127  2.62e-13 ***
#> ---
#> Expl. vars.: poly(Temp, 2)
#>                Estimate Std. Error  z value  Pr(>|z|)
#> (Intercept)      5.1026   0.052341  97.4880 < 2.2e-16 ***
#> poly(Temp, 2)1  -2.6300   0.635659  -4.1374   3.5e-05 ***
#> poly(Temp, 2)2  -5.5977   0.432368 -12.9470 < 2.2e-16 ***
#>
#>
#> # SAMPLE: 6
#>
#>
#> ### Dep. var.: Ozone
#>
#> Expl. vars.: poly(Wind, 2)
#>                Estimate Std. Error z value  Pr(>|z|)
#> (Intercept)      3.0710   0.124756 24.6160 < 2.2e-16 ***
#> poly(Wind, 2)1   6.0990   1.922800  3.1720  0.001514 **
#> poly(Wind, 2)2  -2.2844   1.031900 -2.2138  0.026842 *
#> ---
#> Expl. vars.: poly(Temp, 2)
#>                Estimate Std. Error z value  Pr(>|z|)
#> (Intercept)      3.4184   0.080824 42.2950 < 2.2e-16 ***
#> poly(Temp, 2)1   3.2650   1.293500  2.5242  0.011595 *
#> poly(Temp, 2)2   5.6690   1.562600  3.6279  0.000286 ***
#>
#> ### Dep. var.: Solar.R
#>
#> Expl. vars.: poly(Wind, 2)
#>                 Estimate Std. Error  z value  Pr(>|z|)
#> (Intercept)     5.218500   0.014191 367.7400 < 2.2e-16 ***
#> poly(Wind, 2)1  2.220400   0.177688  12.4960 < 2.2e-16 ***
#> poly(Wind, 2)2 -0.624761   0.139873  -4.4666  7.95e-06 ***
#> ---
#> Expl. vars.: poly(Temp, 2)
#>                Estimate Std. Error    z value  Pr(>|z|)
#> (Intercept)    5.200300   0.019697 264.010000 < 2.2e-16 ***
#> poly(Temp, 2)1 3.286100   0.364969   9.003800 < 2.2e-16 ***
#> poly(Temp, 2)2 0.155763   0.322767   0.482586   0.62939
#>
#>
#> # SAMPLE: 7
#>
#>
#> ### Dep. var.: Ozone
#>
#> Expl. vars.: poly(Wind, 2)
#>                 Estimate Std. Error z value  Pr(>|z|)
#> (Intercept)     3.778200   0.043180 87.4980 < 2.2e-16 ***
#> poly(Wind, 2)1 -6.722500   0.703358 -9.5577 < 2.2e-16 ***
#> poly(Wind, 2)2 -0.959048   0.681576 -1.4071  0.159397
#> ---
#> Expl. vars.: poly(Temp, 2)
#>                Estimate Std. Error z value  Pr(>|z|)
#> (Intercept)       2.595   0.190413 13.6280 < 2.2e-16 ***
#> poly(Temp, 2)1   24.895   3.062000  8.1301  4.29e-16 ***
#> poly(Temp, 2)2   -8.578   1.768700 -4.8499  1.24e-06 ***
#>
#> ### Dep. var.: Solar.R
#>
#> Expl. vars.: poly(Wind, 2)
#>                 Estimate Std. Error    z value  Pr(>|z|)
#> (Intercept)     5.360500   0.014178 378.100000 < 2.2e-16 ***
#> poly(Wind, 2)1 -0.682752   0.216716  -3.150400   0.00163 **
#> poly(Wind, 2)2  0.013888   0.277451   0.050054  0.960079
#> ---
#> Expl. vars.: poly(Temp, 2)
#>                Estimate Std. Error z value  Pr(>|z|)
#> (Intercept)      5.1094   0.056186 90.9360 < 2.2e-16 ***
#> poly(Temp, 2)1   4.8281   0.928155  5.2018  1.97e-07 ***
#> poly(Temp, 2)2  -1.1147   0.615860 -1.8101  0.070285 .
#>
#>
#> # SAMPLE: 8
#>
#>
#> ### Dep. var.: Ozone
#>
#> Expl. vars.: poly(Wind, 2)
#>                 Estimate Std. Error   z value  Pr(>|z|)
#> (Intercept)     3.755800   0.041522 90.454000 < 2.2e-16 ***
#> poly(Wind, 2)1 -6.956500   0.697914 -9.967500 < 2.2e-16 ***
#> poly(Wind, 2)2 -0.344133   0.486153 -0.707869  0.479027
#> ---
#> Expl. vars.: poly(Temp, 2)
#>                Estimate Std. Error z value  Pr(>|z|)
#> (Intercept)      3.1434   0.097610 32.2040 < 2.2e-16 ***
#> poly(Temp, 2)1  17.2990   1.746300  9.9064 < 2.2e-16 ***
#> poly(Temp, 2)2  -5.7801   0.896709 -6.4459  1.15e-10 ***
#>
#> ### Dep. var.: Solar.R
#>
#> Expl. vars.: poly(Wind, 2)
#>                 Estimate Std. Error  z value  Pr(>|z|)
#> (Intercept)     5.134500   0.016224 316.4800 < 2.2e-16 ***
#> poly(Wind, 2)1 -0.459093   0.272152  -1.6869  0.091623 .
#> poly(Wind, 2)2  0.683835   0.246613   2.7729  0.005556 **
#> ---
#> Expl. vars.: poly(Temp, 2)
#>                 Estimate Std. Error  z value  Pr(>|z|)
#> (Intercept)     4.926600   0.041382 119.0500 < 2.2e-16 ***
#> poly(Temp, 2)1  4.222000   0.816701   5.1696  2.35e-07 ***
#> poly(Temp, 2)2 -0.794082   0.462440  -1.7172   0.08595 .
#>
#>
#> # SAMPLE: 9
#>
#>
#> ### Dep. var.: Ozone
#>
#> Expl. vars.: poly(Wind, 2)
#>                Estimate Std. Error z value  Pr(>|z|)
#> (Intercept)      3.3642   0.036953 91.0400 < 2.2e-16 ***
#> poly(Wind, 2)1  -4.0441   0.492839 -8.2057  2.29e-16 ***
#> poly(Wind, 2)2   2.9891   0.470324  6.3554  2.08e-10 ***
#> ---
#> Expl. vars.: poly(Temp, 2)
#>                Estimate Std. Error z value  Pr(>|z|)
#> (Intercept)      3.3695   0.038292 87.9930 < 2.2e-16 ***
#> poly(Temp, 2)1   6.1002   0.690949  8.8287 < 2.2e-16 ***
#> poly(Temp, 2)2   2.2548   0.627045  3.5960  0.000323 ***
#>
#> ### Dep. var.: Solar.R
#>
#> Expl. vars.: poly(Wind, 2)
#>                 Estimate Std. Error    z value  Pr(>|z|)
#> (Intercept)     5.123000   0.014639 349.970000 < 2.2e-16 ***
#> poly(Wind, 2)1 -0.564257   0.193672  -2.913500  0.003574 **
#> poly(Wind, 2)2  0.144413   0.224957   0.641955  0.520902
#> ---
#> Expl. vars.: poly(Temp, 2)
#>                 Estimate Std. Error  z value  Pr(>|z|)
#> (Intercept)     5.111400   0.015302 334.0300 < 2.2e-16 ***
#> poly(Temp, 2)1  1.286100   0.251169   5.1203  3.05e-07 ***
#> poly(Temp, 2)2 -0.928769   0.260074  -3.5712  0.000355 ***
# Different way of displaying the results with "compact"
summary(est_split, "compact")
#>    sample     lhs           rhs      (Intercept)   poly(Wind, 2)1
#> 1       5 Ozone   poly(Wind, 2)  3.3*** (0.0424) -6.41*** (0.617)
#> 2       5 Ozone   poly(Temp, 2) 4.06*** (0.109)
#> 3       5 Solar.R poly(Wind, 2) 5.21*** (0.0182)   -0.57* (0.272)
#> 4       5 Solar.R poly(Temp, 2)  5.1*** (0.0523)
#> 5       6 Ozone   poly(Wind, 2) 3.07*** (0.125)     6.1** (1.92)
#> 6       6 Ozone   poly(Temp, 2) 3.42*** (0.0808)
#> 7       6 Solar.R poly(Wind, 2) 5.22*** (0.0142)  2.22*** (0.178)
#> 8       6 Solar.R poly(Temp, 2)  5.2*** (0.0197)
#> 9       7 Ozone   poly(Wind, 2) 3.78*** (0.0432) -6.72*** (0.703)
#> 10      7 Ozone   poly(Temp, 2)  2.6*** (0.19)
#> 11      7 Solar.R poly(Wind, 2) 5.36*** (0.0142) -0.683** (0.217)
#> 12      7 Solar.R poly(Temp, 2) 5.11*** (0.0562)
#> 13      8 Ozone   poly(Wind, 2) 3.76*** (0.0415) -6.96*** (0.698)
#> 14      8 Ozone   poly(Temp, 2) 3.14*** (0.0976)
#> 15      8 Solar.R poly(Wind, 2) 5.13*** (0.0162)  -0.459. (0.272)
#> 16      8 Solar.R poly(Temp, 2) 4.93*** (0.0414)
#> 17      9 Ozone   poly(Wind, 2) 3.36*** (0.037)  -4.04*** (0.493)
#> 18      9 Ozone   poly(Temp, 2) 3.37*** (0.0383)
#> 19      9 Solar.R poly(Wind, 2) 5.12*** (0.0146) -0.564** (0.194)
#> 20      9 Solar.R poly(Temp, 2) 5.11*** (0.0153)
#>       poly(Wind, 2)2   poly(Temp, 2)1    poly(Temp, 2)2
#> 1    3.08*** (0.57)
#> 2                     11.3*** (1.45)      2.51* (1.17)
#> 3   -1.83*** (0.251)
#> 4                    -2.63*** (0.636)   -5.6*** (0.432)
#> 5     -2.28* (1.03)
#> 6                       3.26* (1.29)    5.67*** (1.56)
#> 7  -0.625*** (0.14)
#> 8                     3.29*** (0.365)     0.156 (0.323)
#> 9     -0.959 (0.682)
#> 10                    24.9*** (3.06)   -8.58*** (1.77)
#> 11    0.0139 (0.277)
#> 12                    4.83*** (0.928)    -1.11. (0.616)
#> 13    -0.344 (0.486)
#> 14                    17.3*** (1.75)   -5.78*** (0.897)
#> 15   0.684** (0.247)
#> 16                    4.22*** (0.817)   -0.794. (0.462)
#> 17   2.99*** (0.47)
#> 18                     6.1*** (0.691)   2.25*** (0.627)
#> 19     0.144 (0.225)
#> 20                    1.29*** (0.251) -0.929*** (0.26)
# You can still select which sample/LHS/RHS to display
est_split[sample = 1:2, lhs = 1, rhs = 1]
#> Standard-errors: Standard
#> Dep. var.: Ozone
#> Expl. vars.: poly(Wind, 2)
#> Sample: 5
#>                Estimate Std. Error  z value  Pr(>|z|)
#> (Intercept)      3.3040   0.042415  77.8970 < 2.2e-16 ***
#> poly(Wind, 2)1  -6.4147   0.617300 -10.3910 < 2.2e-16 ***
#> poly(Wind, 2)2   3.0758   0.570194   5.3944  6.88e-08 ***
#> ---
#> Sample: 6
#>                Estimate Std. Error z value  Pr(>|z|)
#> (Intercept)      3.0710   0.124756 24.6160 < 2.2e-16 ***
#> poly(Wind, 2)1   6.0990   1.922800  3.1720  0.001514 **
#> poly(Wind, 2)2  -2.2844   1.031900 -2.2138  0.026842 *