Nothing
R/coefplot.R
In fixest: Fast Fixed-Effects Estimations
Defines functions
plot.fixest
getFixest_coefplot
setFixest_coefplot
gen_iplot
expr_builder
replace_and_make_callable
coefplot_prms
coefplot
Documented in
coefplot
getFixest_coefplot
plot.fixest
setFixest_coefplot
#----------------------------------------------#
# Author: Laurent Berge
# Date creation: Mon Feb 10 09:46:01 2020
# ~: Simple function to display results
# from multiple estimations
#----------------------------------------------#
#' Plots confidence intervals and point estimates
#'
#' This function plots the results of estimations (coefficients and confidence intervals).
#' The function `iplot` restricts the output to variables created with [`i`], either
#' interactions with factors or raw factors.
#'
#' @inheritParams etable
#' @inheritSection etable Arguments keep, drop and order
#'
#' @param ... Estimation results, which can be of the following form:
#' i) an estimation object (obtained for example from [`feols`]),
#' ii) a matrix of coefficients table, iii) a numeric vector of the point
#' estimates -- the latter requiring the extra arguments `se` or `ci_low` and `ci_high`.
#' @param objects A list of `fixest` estimation objects, or `NULL` (default). If provided,
#' the objects in `...` are ignored and the only coefficients reported are the ones in the
#' argument `objects`.
#' @param vcov Versatile argument to specify the VCOV.
#' In general, it is either a character scalar equal to a VCOV type, either a formula of the form:
#' vcov_type ~ variables. The VCOV types implemented are: "iid", "hetero" (or "HC1"),
#' "cluster", "twoway", "NW" (or "newey_west"), "DK" (or "driscoll_kraay"), and "conley".
#' It also accepts object from vcov_cluster, vcov_NW, NW, vcov_DK, DK, vcov_conley and conley.
#' It also accepts covariance matrices computed externally.
#' Finally it accepts functions to compute the covariances.
#' See the vcov documentation in the vignette.
#'
#' You can pass several VCOVs (as above) if you nest them into a list.
#' If the number of VCOVs equals the number of models, eahc VCOV is mapped to the appropriate model.
#' If there is one model and several VCOVs, or if the first element of the list is equal to
#' `"each"` or `"times"`, then the estimations will be replicated and the results
#' for each estimation and each VCOV will be reported.
#' @param se The standard errors of the estimates. It may be missing.
#' @param ci_low If `se` is not provided, the lower bound of the confidence interval.
#' For each estimate.
#' @param ci_high If `se` is not provided, the upper bound of the confidence interval.
#' For each estimate.
#' @param horiz A logical scalar, default is `FALSE`. Whether to display the confidence
#' intervals horizontally instead of vertically.
#' @param x The value of the x-axis. If missing, the names of the argument `estimate`
#' are used.
#' @param x.shift Shifts the confidence intervals bars to the left or right, depending
#' on the value of `x.shift`. Default is 0.
#' @param ci.width The width of the extremities of the confidence intervals. Default is `0.1`.
#' @param ci_level Scalar between 0 and 1: the level of the CI. By default it is equal to 0.95.
#' @param add Default is `FALSE`, if the intervals are to be added to an existing
#' graph. Note that if it is the case, then the argument `x` MUST be numeric.
#' @param xlim Numeric vector of length 2 which gives the limits of the plotting region for
#' the x-axis. The default is `NULL`, which means that it is automatically defined.
#' Use the argument `xlim.add` to simply increase or decrese the default limits.
#' @param ylim Numeric vector of length 2 which gives the limits of the plotting region for
#' the y-axis. The default is `NULL`, which means that it is automatically defined.
#' Use the argument `ylim.add` to simply increase or decrese the default limits.
#' @param pch The patch of the coefficient estimates. Default is 1 (circle).
#' This is an alias to tha argument `pt.pch`.
#' @param pt.pch The patch of the coefficient estimates. Default is 1 (circle).
#' @param cex Numeric, default is 1. Expansion factor for the points
#' @param pt.cex The size of the coefficient estimates. Default is the other argument `cex`.
#' @param col The color of the points and the confidence intervals. Default is 1
#' ("black"). Note that you can set the colors separately for each of them
#' with `pt.col` and `ci.col`.
#' @param pt.col The color of the coefficient estimates. Default is equal to the argument `col`.
#' @param ci.col The color of the confidence intervals. Default is equal to the argument `col`.
#' @param lwd General line with. Default is 1.
#' @param lty The line type of the confidence intervals. Default is 1.
#' This is an alias to the argument `ci.lty`.
#' @param pt.lwd The line width of the coefficient estimates. Default is equal to
#' the other argument `lwd`.
#' @param ci.lwd The line width of the confidence intervals. Default is equal to
#' the other argument `lwd`.
#' @param ci.lty The line type of the confidence intervals. Default is 1.
#' @param grid Logical, default is `TRUE`. Whether a grid should be displayed. You
#' can set the display of the grid with the argument `grid.par`.
#' @param grid.par List. Parameters of the grid. The default values are: `lty =
#' 3` and `col = "gray"`. You can add any graphical parameter that will be passed
#' to [`graphics::abline`]. You also have two additional arguments: use `horiz =
#' FALSE` to disable the horizontal lines, and use `vert = FALSE` to disable the
#' vertical lines. Eg: `grid.par = list(vert = FALSE, col = "red", lwd = 2)`.
#' @param zero Logical scalar, default is `TRUE`. Whether the 0 should be displayed
#' in the limits of the y-axis.
#' Note that you can set how this zero line looks like with the argument `zero.par`.
#' @param zero.par A named list of graphical parameters or a logical scalar.
#' This argument is a list containing the graphical parameters used to draw the zero-line.
#' The default value is `list(col = "black", lwd = 1)` (it's the same if `TRUE`).
#' Set it to `FALSE` to turn off the special emphasis of the zero line.
#' You can add any graphical parameter that will be passed
#' to [`graphics::abline`]. Example: `zero.par = list(col = "darkblue", lwd = 3)`.
#' @param pt.join Logical, default is `FALSE`. If `TRUE`, then the coefficient estimates
#' are joined with a line.
#' @param df.t Integer scalar or `NULL` (default). The degrees of freedom (DoF) to use
#' when computing the confidence intervals with the Student t. By default it
#' tries to capture the DoF from the estimation. To use a Normal law to compute the
#' confidence interval, use `df.t = Inf`.
#' @param pt.join.par List. Parameters of the line joining the coefficients. The
#' default values are: `col = pt.col` and `lwd = lwd`. You can add any graphical
#' parameter that will be passed to [`lines`]. Eg: `pt.join.par = list(lty = 2)`.
#' @param ref Used to add points at `y = 0` (typically to visualize reference points).
#' Either: i) "auto" (default), ii) a character vector of length 1, iii) a list
#' of length 1, iv) a named integer vector of length 1, or v) a numeric vector.
#' By default, in `iplot`, if the argument `ref` has been used in the estimation,
#' these references are automatically added. If ii), ie a character scalar, then
#' that coefficient equal to zero is added as the first coefficient. If a list or
#' a named integer vector of length 1, then the integer gives the position of the
#' reference among the coefficients and the name gives the coefficient name. A non-named
#' numeric value of `ref` only works if the x-axis is also numeric (which can happen
#' in `iplot`).
#' @param ref.line Logical or numeric, default is "auto", whose behavior depends
#' on the situation. It is `TRUE` only if: i) interactions are plotted, ii) the
#' x values are numeric and iii) a reference is found. If `TRUE`, then a vertical
#' line is drawn at the level of the reference value. Otherwise, if numeric a vertical
#' line will be drawn at that specific value.
#' @param ref.line.par List. Parameters of the vertical line on the reference. The
#' default values are: `col = "black"` and `lty = 2`. You can add any graphical
#' parameter that will be passed to [`graphics::abline`].
#' Eg: `ref.line.par = list(lty = 1, lwd = 3)`.
#' @param xlim.add A numeric vector of length 1 or 2. It represents an extension
#' factor of xlim, in percentage. Eg: `xlim.add = c(0, 0.5)` extends `xlim` of 50%
#' on the right. If of length 1, positive values represent the right, and negative
#' values the left (Eg: `xlim.add = -0.5` is equivalent to `xlim.add = c(0.5, 0)`).
#' @param ylim.add A numeric vector of length 1 or 2. It represents an extension
#' factor of ylim, in percentage. Eg: `ylim.add = c(0, 0.5)` extends `ylim` of 50%
#' on the top. If of length 1, positive values represent the top, and negative values
#' the bottom (Eg: `ylim.add = -0.5` is equivalent to `ylim.add = c(0.5, 0)`).
#' @param only.params Logical, default is `FALSE`. If `TRUE` no graphic is displayed,
#' only the values of `x` and `y` used in the plot are returned.
#' @param ... Other arguments to be passed to `summary`, if `object` is an estimation,
#' and/or to the function `plot` or `lines` (if `add = TRUE`).
#' @param sep The distance between two estimates -- only when argument `object`
#' is a list of estimation results.
#' @param as.multiple Logical: default is `FALSE`. Only when `object` is a single
#' estimation result: whether each coefficient should have a different color, line
#' type, etc. By default they all get the same style.
#' @param bg Background color for the plot. By default it is white.
#' @param ci.join Logical default to `FALSE`. Whether to join the extremities of
#' the confidence intervals. If `TRUE`, then you can set the graphical parameters
#' with the argument `ci.join.par`.
#' @param ci.join.par A list of parameters to be passed to [`graphics::lines`].
#' Only used if `ci.join=TRUE`. By default it is equal to `list(lwd = lwd, col = col, lty = 2)`.
#' @param ci.fill Logical default to `FALSE`. Whether to fill the confidence intervals
#' with a color. If `TRUE`, then you can set the graphical parameters
#' with the argument `ci.fill.par`.
#' @param ci.fill.par A list of parameters to be passed to [`graphics::polygon`].
#' Only used if `ci.fill=TRUE`. By default it is equal to `list(col = "lightgray", alpha = 0.5)`.
#' Note that `alpha` is a special parameter that adds transparency to the color
#' (ranges from 0 to 1).
#' @param group A list, default is missing. Each element of the list reports the
#' coefficients to be grouped while the name of the element is the group name. Each
#' element of the list can be either: i) a character vector of length 1, ii) of
#' length 2, or ii) a numeric vector. If equal to: i) then it is interpreted as
#' a pattern: all element fitting the regular expression will be grouped (note that
#' you can use the special character "^^" to clean the beginning of the names, see
#' example), if ii) it corresponds to the first and last elements to be grouped,
#' if iii) it corresponds to the coefficients numbers to be grouped. If equal to
#' a character vector, you can use a percentage to tell the algorithm to look at
#' the coefficients before aliasing (e.g. `"%varname"`). Example of valid uses:
#' `group=list(group_name=\"pattern\")`, `group=list(group_name=c(\"var_start\", \"var_end\"))`,
#' `group=list(group_name=1:2))`. See details.
#' @param group.par A list of parameters controlling the display of the group. The
#' parameters controlling the line are: `lwd`, `tcl` (length of the tick), `line.adj`
#' (adjustment of the position, default is 0), `tick` (whether to add the ticks),
#' `lwd.ticks`, `col.ticks`. Then the parameters controlling the text: `text.adj`
#' (adjustment of the position, default is 0), `text.cex`, `text.font`, `text.col`.
#' @param pt.bg The background color of the point estimate (when the `pt.pch` is
#' in 21 to 25). Defaults to NULL.
#' @param lab.cex The size of the labels of the coefficients. Default is missing.
#' It is automatically set by an internal algorithm which can go as low as `lab.min.cex`
#' (another argument).
#' @param lab.min.cex The minimum size of the coefficients labels, as set by the
#' internal algorithm. Default is 0.85.
#' @param lab.max.mar The maximum size the left margin can take when trying to fit
#' the coefficient labels into it (only when `horiz = TRUE`). This is used in the
#' internal algorithm fitting the coefficient labels. Default is `0.25`.
#' @param lab.fit The method to fit the coefficient labels into the plotting region
#' (only when `horiz = FALSE`). Can be `"auto"` (the default), `"simple"`, `"multi"`
#' or `"tilted"`. If `"simple"`, then the classic axis is drawn. If `"multi"`, then
#' the coefficient labels are fit horizontally across several lines, such that they
#' don't collide. If `"tilted"`, then the labels are tilted. If `"auto"`, an automatic
#' choice between the three is made.
#' @param main The title of the plot. Default is `"Effect on __depvar__"`. You can
#' use the special variable `__depvar__` to set the title (useful when you set the
#' plot default with [`setFixest_coefplot`]).
#' @param value.lab The label to appear on the side of the coefficient values. If
#' `horiz = FALSE`, the label appears in the y-axis. If `horiz = TRUE`, then it
#' appears on the x-axis. The default is equal to `"Estimate and __ci__ Conf. Int."`,
#' with `__ci__` a special variable giving the value of the confidence interval.
#' @param xlab The label of the x-axis, default is `NULL`. Note that if `horiz =
#' TRUE`, it overrides the value of the argument `value.lab`.
#' @param ylab The label of the y-axis, default is `NULL`. Note that if
#' `horiz = FALSE`, it overrides the value of the argument `value.lab`.
#' @param sub A subtitle, default is `NULL`.
#' @param style A character scalar giving the style of the plot to be used. You
#' can set styles with the function [`setFixest_coefplot`], setting all the default
#' values of the function. If missing, then it switches to either "default" or "iplot",
#' depending on the calling function.
#' @param i.select Integer scalar, default is 1. In `iplot`, used to select which
#' variable created with `i()` to select. Only used when there are several variables
#' created with [`i`]. This is an index, just try increasing numbers to hopefully
#' obtain what you want. Note that it works much better when the variables are "pure"
#' `i()` and not interacted with other variables. For example: `i(species, x1)`
#' is good while `i(species):x1` isn't. The latter will also work but the index
#' may feel weird in case there are many `i()` variables.
#' @param do_iplot Logical, default is `FALSE`. For internal use only.
#' If `TRUE`, then `iplot` is run instead of `coefplot`.
#' @param plot_prms A named list. It may contain additionnal parameters to be passed
#' to the plot.
#'
#' @seealso
#' See [`setFixest_coefplot`] to set the default values of `coefplot`, and the estimation
#' functions: e.g. [`feols`], [`fepois`][fixest::feglm], [`feglm`], [`fenegbin`][fixest::femlm].
#'
#' @section Setting custom default values:
#' The function `coefplot` dispose of many arguments to parametrize the plots. Most
#' of these arguments can be set once an for all using the function [`setFixest_coefplot`].
#' See Example 3 below for a demonstration.
#'
#' @section iplot:
#'
#' The function `iplot` restricts `coefplot` to interactions or factors created
#' with the function [`i`]. Only *one* of the i-variables will be plotted at a time.
#' If you have several i-variables, you can navigate through them with the `i.select` argument.
#'
#' The argument `i.select` is an index that will go through all the i-variables.
#' It will work well if the variables are pure, meaning not interacted with other
#' variables. If the i-variables are interacted, the index may have an odd behavior
#' but will (in most cases) work all the same, just try some numbers up until you
#' (hopefully) obtain the graph you want.
#'
#' Note, importantly, that interactions of two factor variables are (in general)
#' disregarded since they would require a 3-D plot to be properly represented.
#'
#' @section Mathematical expressions:
#'
#' You can add [`plotmath`][grDevices::plotmath] mathematical expressions in the arguments
#' `main`, `sub`, `xlab`, or `ylab`. To do so, start the character string with an ampersand.
#' For example main = "&lambda^2".
#'
#'
#' @author
#' Laurent Berge
#'
#' @examples
#'
#' #
#' # Example 1: Stacking two sets of results on the same graph
#' #
#'
#' # Estimation on Iris data with one fixed-effect (Species)
#' # + we cluster the standard-errors
#' est = feols(Petal.Length ~ Petal.Width + Sepal.Width | Species,
#' iris, vcov = "cluster")
#'
#' # Now with "regular" standard-errors
#' est_std = summary(est, vcov = "iid")
#'
#' # You can plot the two results at once
#' coefplot(est, est_std)
#'
#' # You could also use the argument vcov
#' coefplot(est, vcov = list("cluster", "iid"))
#'
#' # Alternatively, you can use the argument x.shift
#' # to do it sequentially:
#'
#' # First graph with clustered standard-errors
#' coefplot(est, x.shift = -.2)
#'
#' # 'x.shift' was used to shift the coefficients to the left.
#'
#' # Second set of results: this time with
#' # standard-errors that are not clustered.
#' coefplot(est, vcov = "iid", x.shift = .2,
#' add = TRUE, col = 2, ci.lty = 2, pch = 15)
#'
#' legend("topright", col = 1:2, pch = 20, lwd = 1, lty = 1:2,
#' legend = c("Clustered", "IID"), title = "Standard-Errors")
#'
#'
#' #
#' # Example 2: Interactions
#' #
#'
#'
#' # Now we estimate and plot the "yearly" treatment effects
#'
#' data(base_did)
#' base_inter = base_did
#'
#' # We interact the variable 'period' with the variable 'treat'
#' est_did = feols(y ~ x1 + i(period, treat, 5) | id + period, base_inter)
#'
#' # In the estimation, the variable treat is interacted
#' # with each value of period but 5, set as a reference
#'
#' # coefplot will show all the coefficients:
#' coefplot(est_did)
#'
#' # Note that the grouping of the coefficients is due to 'group = "auto"'
#'
#' # If you want to keep only the coefficients
#' # created with i() (ie the interactions), use iplot
#' iplot(est_did)
#'
#' # We can see that the graph is different from before:
#' # - only interactions are shown,
#' # - the reference is present,
#' # => this is fully flexible
#'
#' iplot(est_did, ref.line = FALSE, pt.join = TRUE)
#'
#'
#' #
#' # What if the interacted variable is not numeric?
#'
#' # Let's create a "month" variable
#' all_months = c("aug", "sept", "oct", "nov", "dec", "jan",
#' "feb", "mar", "apr", "may", "jun", "jul")
#' base_inter$period_month = all_months[base_inter$period]
#'
#' # The new estimation
#' est = feols(y ~ x1 + i(period_month, treat, "oct") | id+period, base_inter)
#' # Since 'period_month' of type character, coefplot sorts it
#' iplot(est)
#'
#' # To respect a plotting order, use a factor
#' base_inter$month_factor = factor(base_inter$period_month, levels = all_months)
#' est = feols(y ~ x1 + i(month_factor, treat, "oct") | id + period, base_inter)
#' iplot(est)
#'
#'
#' #
#' # Example 3: Setting defaults
#' #
#'
#' # coefplot has many arguments, which makes it highly flexible.
#' # If you don't like the default style of coefplot. No worries,
#' # you can set *your* default by using the function
#' # setFixest_coefplot()
#'
#' dict = c("Petal.Length"="Length (Petal)", "Petal.Width"="Width (Petal)",
#' "Sepal.Length"="Length (Sepal)", "Sepal.Width"="Width (Sepal)")
#'
#' setFixest_coefplot(ci.col = 2, pt.col = "darkblue", ci.lwd = 3,
#' pt.cex = 2, pt.pch = 15, ci.width = 0, dict = dict)
#'
#' est = feols(Petal.Length ~ Petal.Width + Sepal.Length +
#' Sepal.Width + i(Species), iris)
#'
#' # And that's it
#' coefplot(est)
#'
#' # You can set separate default values for iplot
#' setFixest_coefplot("iplot", pt.join = TRUE, pt.join.par = list(lwd = 2, lty = 2))
#' iplot(est)
#'
#' # To reset to the default settings:
#' setFixest_coefplot("all", reset = TRUE)
#' coefplot(est)
#'
#' #
#' # Example 4: group + cleaning
#' #
#'
#' # You can use the argument group to group variables
#' # You can further use the special character "^^" to clean
#' # the beginning of the coef. name: particularly useful for factors
#'
#' est = feols(Petal.Length ~ Petal.Width + Sepal.Length +
#' Sepal.Width + Species, iris)
#'
#' # No grouping:
#' coefplot(est)
#'
#' # now we group by Sepal and Species
#' coefplot(est, group = list(Sepal = "Sepal", Species = "Species"))
#'
#' # now we group + clean the beginning of the names using the special character ^^
#' coefplot(est, group = list(Sepal = "^^Sepal.", Species = "^^Species"))
#'
#'
coefplot = function(..., objects = NULL, style = NULL, se, ci_low, ci_high, df.t = NULL,
vcov = NULL, cluster = NULL,
x, x.shift = 0, horiz = FALSE,
dict = NULL, keep, drop, order, ci.width = "1%",
ci_level = 0.95, add = FALSE, plot_prms = list(),
pch = c(20, 17, 15, 21, 24, 22), col = 1:8, cex = 1, lty = 1, lwd = 1,
ylim = NULL, xlim = NULL,
pt.pch = pch,
pt.bg = NULL, pt.cex = cex, pt.col = col, ci.col = col, pt.lwd = lwd,
ci.lwd = lwd, ci.lty = lty, grid = TRUE,
grid.par = list(lty = 3, col = "gray"),
zero = TRUE, zero.par = list(col = "black", lwd = 1), pt.join = FALSE,
pt.join.par = list(col = pt.col, lwd = lwd), ci.join = FALSE,
ci.join.par = list(lwd = lwd, col = col, lty = 2), ci.fill = FALSE,
ci.fill.par = list(col = "lightgray", alpha = 0.5), ref = "auto",
ref.line = "auto", ref.line.par = list(col = "black", lty = 2), lab.cex,
lab.min.cex = 0.85, lab.max.mar = 0.25, lab.fit = "auto", xlim.add,
ylim.add, only.params = FALSE, sep, as.multiple = FALSE,
bg, group = "auto", group.par = list(lwd = 2, line = 3, tcl = 0.75),
main = "Effect on __depvar__", value.lab = "Estimate and __ci__ Conf. Int.",
ylab = NULL, xlab = NULL, sub = NULL, i.select = NULL, do_iplot = NULL){
check_set_arg(lab.fit, "match(auto, simple, multi, tilted)")
dots = list(...)
ylab_add_ci = missing(ci_low)
#
# iplot
#
# deprecated sd argument
if("sd" %in% names(dots)){
se = dots$sd
dots$sd = NULL
}
# catching the deprecated `internal.only.i` argument
if("internal.only.i" %in% names(dots)){
is_iplot = isTRUE(dots$internal.only.i)
dots$internal.only.i = NULL
} else {
is_iplot = isTRUE(do_iplot)
}
if(is_iplot){
check_arg(i.select, "NULL integer scalar GE{1}")
if(is.null(i.select)){
i.select = 1
}
}
if(!is.null(objects)){
if(!is.list(objects) || inherits(objects, "fixest")){
dots = list(objects)
} else {
dots = objects
}
}
if(length(dots) == 0){
stop_up("You must provide at least one model to plot. Currently no model is provided.",
up = 1 * is_iplot)
}
# retro-compatibility
if("object" %in% names(dots)){
all_dots = dots$object
dots_rest = dots[names(dots) != "object"]
if(!is.null(oldClass(all_dots))){
all_dots = append(list(all_dots), dots_rest)
} else {
all_dots = append(all_dots, dots_rest)
}
} else {
all_dots = dots
}
info_models = flatten_list_of_models(all_dots, accept_non_fixest = TRUE)
all_models = info_models$all_models
#
# multiple VOCVs
#
# we only keep the arg `vcov`
if(!missnull(cluster)){
vcov = oldargs_to_vcov(NULL, cluster, vcov)
}
if(is.list(vcov) && length(vcov) > 1){
# we multiply the models
is_fixest = sapply(all_models, inherits, "fixest")
if(!all(is_fixest)){
stop_up("To use the `vcov` as a list, and hence have different VCOVs ",
"for different models, you provide only `fixest` estimations.",
"\nProblem: the {nth ? which(!is_fixest)} objects are not `fixest` models.")
}
vcov_1 = vcov[[1]]
is_rep = identical(vcov_1, "times") || identical(vcov_1, "each")
if(length(all_models) > 1 && !is_rep){
stop_up("If 'vcov' is a list, it must be of the same length as the number of models, ",
"or you should add the 'each' or 'times' keyword as the first ",
"element of the list.",
"\nProblem: there are {len ? all_models} models vs {len ? vcov} elements in `vcov`.")
}
n_models = length(all_models)
all_vcov = vcov
if(is_rep || n_models == 1){
if(is_rep){
all_vcov[[1]] = NULL
}
n_vcov = length(all_vcov)
n_total = n_models * n_vcov
if(vcov_1 == "times"){
id_mod = rep(1:n_models, n_vcov)
id_vcov = rep(1:n_vcov, each = n_models)
} else {
IS_EACH = TRUE
id_mod = rep(1:n_models, each = n_vcov)
id_vcov = rep(1:n_vcov, n_models)
}
mega_models = vector("list", n_total)
for(i in 1:n_total){
mega_models[[i]] = summary(all_models[[id_mod[i]]], vcov = all_vcov[[id_vcov[i]]])
}
all_models = mega_models
} else {
for(i in seq_along(vcov)){
all_models[[i]] = summary(all_models[[i]] , vcov = all_vcov[[i]])
}
}
vcov = NULL
}
#
# Setting the default values ####
#
opts = getOption("fixest_coefplot")
check_arg(style, "NULL character scalar")
if(is.null(style)){
if(is_iplot){
style = "iplot"
} else {
style = "default"
}
}
old_style = style
style = try(match.arg(style, choices = names(opts)), silent = TRUE)
if(inherits(style, "try-error")){
warning("Unknow style '", old_style, "', using default style instead.")
style = "default"
}
if(style != "default" && !is.null(opts[["default"]])){
# there is always inheritance from the default style
my_opt = opts[["default"]]
opts_style = opts[[style]]
if(!is.null(opts_style)){
my_opt[names(opts_style)] = opts_style
}
} else {
my_opt = opts[[style]]
}
if(!is.list(my_opt)){
warning("The default values of coefplot for style '", style, "' are ill-formed and therefore reset. Use only setFixest_coefplot for setting the default values.")
setFixest_coefplot(style = style, reset = TRUE)
} else if(length(my_opt) > 0){
arg_2_add = setdiff(names(opts[["default"]]), names(my_opt))
if(length(arg_2_add) > 0){
my_opt[arg_2_add] = opts[["default"]][arg_2_add]
# we reorder
arg_list = names(formals(coefplot))
my_fact = factor(names(my_opt), levels = arg_list)
my_opt = my_opt[base::order(my_fact)]
}
if(is_iplot){
sysOrigin = sys.parent()
mc = match.call(definition = sys.function(sysOrigin), call = sys.call(sysOrigin))
} else {
mc = match.call()
}
arg2set = setdiff(names(my_opt), names(mc))
for(arg in arg2set){
my_arg = my_opt[[arg]]
if(is.name(my_arg) || is.call(my_arg)){
if(length(my_opt$extra_values) > 0){
assign(arg, eval(my_arg, my_opt$extra_values))
} else {
assign(arg, eval(my_arg))
}
} else {
assign(arg, my_arg)
}
}
}
# Set up the dictionary
dict = setup_dict(dict, check = TRUE)
if(!is.null(dict) && any(grepl("^&", names(dict)))){
# Speficic markup to identify coefplot aliases
dict_amp = dict[grepl("^&", names(dict))]
names(dict_amp) = gsub("^&", "", names(dict_amp))
dict[names(dict_amp)] = dict_amp
}
#
# Getting the parameters => function coefplot_prms
#
info = coefplot_prms(all_models = all_models, is_root = TRUE, vcov = vcov,
se = se, ci_low = ci_low, ci_high = ci_high,
x = x, x.shift = x.shift, dict = dict, keep = keep, drop = drop,
order = order, ci_level = ci_level, df.t = df.t, ref = ref,
i.select = i.select,
is_iplot = is_iplot, sep = sep, as.multiple = as.multiple)
prms = info$prms
AXIS_AS_NUM = info$num_axis
x_at = info$at
x_labels = info$labels
x_labels_raw = info$x_labels_raw
varlist = info$varlist
my_xlim = info$xlim
suggest_ref_line = info$suggest_ref_line
multiple_est = info$multiple_est
if(only.params){
return(list(prms = prms, is_iplot = is_iplot, at = x_at, labels = x_labels))
}
check_arg(plot_prms, "named list")
ci_low = prms$ci_low
ci_high = prms$ci_high
x_value = prms$x
if(horiz){
# We just reverse the x/y
tmp = prms$x
prms$x = prms$y
prms$y = tmp
}
check_arg(xlab, "NULL character vector")
if(is.null(xlab) && is_iplot){
xlab = "__i__"
}
#
# Title ####
#
check_arg(xlim, ylim, "NULL numeric vector len(2)")
# xlab / main / ylab / sub
check_arg(value.lab, main, "character scalar")
check_arg(xlab, ylab, sub, "character scalar NULL")
if(horiz){
if(is.null(xlab)) xlab = value.lab
} else {
if(is.null(ylab)) ylab = value.lab
}
if(is.null(xlab)) xlab = ""
if(is.null(ylab)) ylab = ""
if(is.null(sub)) sub = ""
main = replace_and_make_callable(main, varlist)
ylab = replace_and_make_callable(ylab, varlist)
xlab = replace_and_make_callable(xlab, varlist)
sub = replace_and_make_callable(sub, varlist)
main = expr_builder(main)
ylab = expr_builder(ylab)
xlab = expr_builder(xlab)
sub = expr_builder(sub)
plot_prms$main = ""
plot_prms$ylab = ""
plot_prms$xlab = ""
plot_prms$sub = ""
#
# group = auto + Height ####
#
# The value of group = "auto" => renaming the labels
if(identical(group, "auto") && is_iplot == FALSE){
# we change the names of interactions
qui = grepl(":", x_labels_raw)
if(any(qui)){
# we only keep the ones with ':', or '::' or '::' and ':', no more
n_colons = lengths(strsplit(x_labels_raw, ':'))
qui = qui & n_colons <= 3
all_vars_to_group = character()
if(any(qui)){
x_inter = gsub(":.+", "", x_labels_raw[qui])
tx_inter = table(x_inter)
all_vars_to_group = names(tx_inter)[tx_inter >= 2]
}
group = list()
for(i in seq_along(all_vars_to_group)){
var_left = all_vars_to_group[i]
qui_select = substr(x_labels_raw, 1, nchar(var_left)) == var_left
# we require to be next to each other
if(!all(diff(which(qui_select)) == 1)) next
x_select = x_labels_raw[qui_select]
is_inter = TRUE
n_trim = 0
group_regex = NULL
if(all(grepl("[^:]:[^:].*::", x_select))){
# treat:period::1
first_part = strsplit(x_select[1], "::")[[1]][1]
var_right = gsub(".+:", "", first_part)
n_max = nchar(first_part) + 2
} else if(all(grepl("::.*[^:]:[^:]", x_select))){
# period::1:treat
# Note that we always put 'treat' on the left
second_part = strsplit(x_select[1], "::")[[1]][2]
var_right = var_left
var_left = gsub(".+:", "", second_part)
n_max = nchar(var_right) + 2
n_trim = nchar(var_left) + 1
group_regex = paste0("%", escape_regex(var_right), "::.+:", escape_regex(var_left))
} else if(all(grepl("::", x_select))) {
is_inter = FALSE
# This is a fixest factor
n_max = nchar(var_left) + 2
} else {
# we need to find out by ourselves...
n = min(nchar(x_select[1:2]))
x_split = strsplit(substr(x_select[1:2], 1, n), "")
start_diff = which.max(x_split[[1]] != x_split[[2]])
ok = FALSE
for(n_max in n:(nchar(var_left) + 2)){
if(all(grepl(substr(x_select[1], 1, n_max), x_select, fixed = TRUE))){
ok = TRUE
break
}
}
if(!ok) next
var_right = gsub(".+:", "", substr(x_select[1], 1, n_max))
}
if(is_inter){
v_name = dict_apply(c(var_left, var_right), dict)
group_name = replace_and_make_callable("__x__ %*% (__y__ == ldots)",
list(x = v_name[1], y = v_name[2]),
text_as_expr = TRUE)
if(is.null(group_regex)){
group[[group_name]] = escape_regex(paste0("%", var_left, ":", var_right))
} else {
group[[group_name]] = group_regex
}
} else {
v_name = dict_apply(var_left, dict)
group_name = replace_and_make_callable("__x__", list(x = v_name),
text_as_expr = TRUE)
group[[group_name]] = paste0("%^", escape_regex(var_left))
}
# We update the labels
x_labels[qui_select] = substr(x_select, n_max + 1, nchar(x_select) - n_trim)
}
}
}
IS_GROUP = !identical(group, "auto") && !add && !missing(group) && !is.null(group) && length(group) > 0 && !is.null(x_labels)
line_height = par("mai")[1] / par("mar")[1]
if(IS_GROUP){
# This means there are groups!
# We compute the height of the group based on the parameters
# passed in group.par
# => it will be used to adjust the margins
# we need: tcl + text.adj
tcl = 0.75
if(!is.null(group.par$tcl)){
tcl = max(group.par$tcl)
}
line.adj = 0
if(!is.null(group.par$line.adj)){
line.adj = max(group.par$line.adj)
}
text.adj = 0
if(!is.null(group.par$text.adj)){
text.adj = max(group.par$text.adj)
}
text.cex = 1
if(!is.null(group.par$text.cex)){
text.cex = max(group.par$text.cex)
}
group.height = (0.5 + tcl + text.cex + text.adj + line.adj)
# We perform basic checks on the group
if(!is.list(group)) stop("Argument 'group' must be a list.")
# We just take care of the special group + cleaning case
# when group starts with ^^
for(i in seq_along(group)){
my_group = group[[i]]
if(! (length(my_group) == 1 && is.character(my_group) && substr(my_group, 1, 2) == "^^") ) next
if(grepl("^%", my_group)){
qui = grepl(gsub("^%", "", my_group), x_labels_raw)
} else {
qui = grepl(my_group, x_labels)
}
if(!any(qui)){
warning("In argument 'group', the pattern: \"", my_group, "\", did not match any coefficient name.")
group[[i]] = NULL
if(length(group) == 0) IS_GROUP = FALSE
next
}
group[[i]] = range(which(qui))
x_labels = gsub(my_group, "", x_labels)
}
} else {
group.height = 0
}
#
# The limits ####
#
# xlim
if(!missnull(xlim)){
check_arg(xlim, "numeric vector len(2) no na")
my_xlim = xlim
}
if(!missnull(xlim.add)){
if((!is.numeric(xlim.add) || !length(xlim.add) %in% 1:2)){
stop("Argument 'xlim.add' must be a numeric vector of length 1 or 2. It represents an extension factor of xlim, in percentage. (Eg: xlim.add = c(0, 0.5) extends xlim of 50% on the right.) If of lentgh 1, positive values represent the right, and negative values the left (Eg: xlim.add = -0.5 is equivalent to xlim.add = c(0.5, 0)).")
}
if(length(xlim.add) == 1){
if(xlim.add > 0) {
xlim.add = c(0, xlim.add)
} else {
xlim.add = c(xlim.add, 0)
}
}
x_width = diff(my_xlim)
my_xlim = my_xlim + xlim.add * x_width
}
# ylim
my_ylim = range(c(ci_low, ci_high))
if(!missnull(ylim)){
check_arg(ylim, "numeric vector len(2) no na")
my_ylim = ylim
} else if(isTRUE(zero)){
# we add 0 to xlim
if(0 > my_ylim[2]){
my_ylim[2] = 0
} else if(0 < my_ylim[1]){
my_ylim[1] = 0
}
}
if(!missnull(ylim.add)){
if((!length(ylim.add) %in% 1:2 || !is.numeric(ylim.add))){
stop("Argument 'ylim.add' must be a numeric vector of length 1 or 2. It represents an extension factor of ylim, in percentage. (Eg: ylim.add = c(0, 0.5) extends ylim of 50% on the top.) If of lentgh 1, positive values represent the top, and negative values the bottom (Eg: ylim.add = -0.5 is equivalent to ylim.add = c(0.5, 0)).")
}
if(length(ylim.add) == 1){
if(ylim.add > 0) {
ylim.add = c(0, ylim.add)
} else {
ylim.add = c(ylim.add, 0)
}
}
y_width = diff(my_ylim)
my_ylim = my_ylim + ylim.add * y_width
}
#
# Margins and cex setting ####
#
if(!missing(lab.cex)){
lab.min.cex = lab.cex
} else {
lab.cex = 1
}
if(horiz){
# we make the labels fit into the margin
xlab.line = 3
sub.line = 4
# Algorithm:
# - if lab.cex is provided:
# => we fit the labels into the margin // we extend the margin until it fits
# - if lab.cex is NOT provided:
# => we extend the margin so that the label fit. If the margin exceeds lab.max.mar% of the plot space, we reduce the cex
# etc, until it fits. If it still does not fit => we extend the margin until it fits.
#
nlines = group.height + 2
# The last 2 means 1 line on each side of the label
if(lab.cex > lab.min.cex){
# No algorithm
lab.width_in = max(strwidth(x_labels, units = "in", cex = lab.cex))
} else {
# Algorithm
max_mar_width_in = par("pin")[1] * lab.max.mar
while(nlines * line_height + lab.width_in > max_mar_width_in && lab.cex > lab.min.cex){
lab.cex = 0.95 * lab.cex
lab.width_in = max(strwidth(x_labels, units = "in", cex = lab.cex))
}
# Final step => if we go too far, we etend the margin, even beyond max_mar_width_in
if(lab.cex < lab.min.cex){
lab.cex = lab.min.cex
lab.width_in = max(strwidth(x_labels, units = "in", cex = lab.cex))
}
}
group.baseline = 2 + lab.width_in / line_height
nlines = nlines + lab.width_in / line_height
ylab.line = nlines
if(ylab != ""){
nlines = nlines + ceiling(strheight(ylab, units = "in") / line_height)
}
total_width = nlines * line_height
if(total_width > par("mai")[2]){
new_mai = par("mai")
new_mai[2] = total_width + 0.05
op = par(mai = new_mai)
on.exit(par(op))
}
my_xlim = rev(my_xlim)
} else {
# We adjust the margin only if there are groups and they
# don't fit in the original margin
# or if there is a xlab and groups at the same time
ylab.line = 3
LINE_MIN_TILTED = 0.25
LINE_MAX_TILTED = 3
if(lab.fit == "auto"){
# We want to display ALL labels
in_to_usr = diff(my_xlim) / par("pin")[1]
w_all = strwidth(x_labels, cex = lab.cex, units = "in") * in_to_usr
em = strwidth("M", units = "in") * in_to_usr
is_collided = ((w_all[-1] + w_all[-length(w_all)]) / 2 + em) > 1
while(any(is_collided) && lab.cex > lab.min.cex){
lab.cex = 0.95 * lab.cex
w_all = strwidth(x_labels, cex = lab.cex, units = "in") * in_to_usr
is_collided = ((w_all[-1] + w_all[-length(w_all)]) / 2 + em) > 1
}
# switch to multi lines
if(any(is_collided) || lab.cex < lab.min.cex){
lab.info = xaxis_labels(at = x_at, labels = x_labels, line.max = 1,
minCex = lab.min.cex, trunc = Inf, only.params = TRUE,
xlim = my_xlim)
if(lab.info$failed){
# switch to tilted
lab.fit = "tilted"
lab.info = xaxis_biased(at = x_at, labels = x_labels,
line.min = LINE_MIN_TILTED,
line.max = LINE_MAX_TILTED,
cex = seq(lab.cex, lab.min.cex, length.out = 4),
trunc = Inf, only.params = TRUE,
ylim = my_ylim)
lab.cex = lab.info$cex
nlines = lab.info$height_line + LINE_MIN_TILTED
} else {
lab.fit = "multi"
nlines = lab.info$height_line
lab.cex = lab.info$cex
}
} else {
lab.fit = "simple"
nlines = 2 * lab.cex
}
} else if(lab.fit == "multi"){
lab.info = xaxis_labels(at = x_at, labels = x_labels, line.max = 1,
minCex = lab.min.cex, trunc = Inf,
only.params = TRUE, xlim = my_xlim)
lab.cex = lab.info$cex
nlines = lab.info$height_line
if(length(unique(lab.info$line)) == 1){
lab.fit = "simple"
nlines = 2 * lab.cex
}
} else if(lab.fit == "tilted"){
lab.info = xaxis_biased(at = x_at, labels = x_labels, line.min = LINE_MIN_TILTED,
line.max = LINE_MAX_TILTED,
cex = seq(lab.cex, lab.min.cex, length.out = 4),
trunc = Inf, only.params = TRUE, ylim = my_ylim)
lab.cex = lab.info$cex
nlines = lab.info$height_line + LINE_MIN_TILTED
} else if(lab.fit == "simple"){
nlines = 2 * lab.cex
}
if(IS_GROUP){
# tcl was set before
group.baseline = nlines + 1 - 0.75 + tcl
nlines = nlines + group.height
xlab.line = nlines
} else {
xlab.line = 3
}
if(xlab != ""){
xlab.line = xlab.line + ceiling(strheight(xlab, units = "in") / line_height) - 1
}
sub.line = xlab.line + 1
if(sub != ""){
nlines = sub.line + 1
} else if(xlab != ""){
nlines = sub.line
}
total_height = nlines * line_height
if(total_height > par("mai")[1]){
new_mai = par("mai")
new_mai[1] = total_height + 0.05
op = par(mai = new_mai)
on.exit(par(op))
}
}
#
# Plot (start) ####
#
all_plot_args = unique(c(names(par()), names(formals(plot.default))))
pblm = setdiff(names(plot_prms), all_plot_args)
if(length(pblm) > 0){
plot_prms[pblm] = NULL
}
# preparation of the do.call
plot_prms$col = col
if(horiz){
plot_prms$xlim = my_ylim
plot_prms$ylim = my_xlim
} else {
plot_prms$xlim = my_xlim
plot_prms$ylim = my_ylim
}
plot_prms$x = prms$x
plot_prms$y = prms$y
plot_prms$type = "p"
#
# Plot Call ####
#
if(!add){
plot_prms$axes = FALSE
# Nude graph
first.par = plot_prms
first.par$type = "n"
do.call("plot", first.par)
# The background
if(!missing(bg)){
dx = diff(my_xlim)
dy = diff(my_ylim)
if(horiz){
rect(xleft = my_ylim[1] - dy, ybottom = my_xlim[1] - dx, xright = my_ylim[2] + dy,
ytop = my_xlim[2] + dx, col = bg)
} else {
rect(xleft = my_xlim[1] - dx, ybottom = my_ylim[1] - dy, xright = my_xlim[2] + dx,
ytop = my_ylim[2] + dy, col = bg)
}
}
# The grid
if(grid){
listDefault(grid.par, "col", "gray")
listDefault(grid.par, "lty", 3)
listDefault(grid.par, "vert", TRUE)
listDefault(grid.par, "horiz", TRUE)
vert = grid.par$vert
g.horiz = grid.par$horiz
grid.par$vert = grid.par$horiz = NULL
if(g.horiz){
do.call("hgrid", grid.par)
}
if(vert){
do.call("vgrid", grid.par)
}
}
check_arg(zero.par, "logical scalar | named list")
if(!isFALSE(zero.par)){
if(isTRUE(zero.par)){
zero.par = list(col = "black", lwd = 1)
}
listDefault(zero.par, "lwd", 1)
listDefault(zero.par, "col", "black")
if(horiz){
zero.par$v = 0
} else {
zero.par$h = 0
}
do.call("abline", zero.par)
}
# Tha annotations
if(main != ""){
title(main = main)
}
if(xlab != ""){
title(xlab = xlab, line = xlab.line)
}
if(ylab != ""){
title(ylab = ylab, line = ylab.line)
}
if(sub != ""){
title(sub = sub, line = sub.line)
}
# Reference line
check_arg(ref.line, "logical scalar | numeric vector no na | charin(auto)")
if(identical(ref.line, "auto")){
ref.line = suggest_ref_line && length(unique(prms[prms$is_ref, "estimate_names_raw"])) == 1
}
is_ref_line = !isFALSE(ref.line)
line_direction = if(horiz) "h" else "v"
if(is.logical(ref.line)){
if(ref.line) {
ref_pblm = is.null(prms$is_ref) || !any(prms$is_ref)
if(ref_pblm && !"v" %in% names(ref.line.par)){
warning("You can use the argument 'ref.line = TRUE' only when a 'natural' reference is found, or if you provided a reference with argument 'ref'. You can still draw vertical lines by using 'v' in argument 'ref.line.par'. Example: ref.line.par=list(v = ", round(x_value[floor(length(x_value)/2)]), ", col=2).")
} else if(!ref_pblm){
where = tapply(prms[prms$is_ref, "x"], prms[prms$is_ref, "estimate_names_raw"], mean)
listDefault(ref.line.par, line_direction, where)
}
}
} else {
listDefault(ref.line.par, line_direction, ref.line)
}
if(is_ref_line){
listDefault(ref.line.par, "lty", 2)
do.call("abline", ref.line.par)
}
box(bty = par("bty"))
if(horiz){
axis(1, lwd = 0, lwd.ticks = 1)
if(AXIS_AS_NUM){
axis(2, las = 1, lwd = 0, lwd.ticks = 1)
} else {
if(any(grepl("^&", x_labels))){
# means we call expression()
# drawback => expressions can overlap
qui = grepl("^&", x_labels)
if(any(!qui)){
axis(2, at = x_at[!qui], labels = x_labels[!qui], las = 1, cex = lab.cex,
lwd = 0, lwd.ticks = 1)
}
for(i in which(qui)){
axis(2, at = x_at[i], labels = expr_builder(x_labels[i]), las = 1, cex = lab.cex,
lwd = 0, lwd.ticks = 1)
}
} else {
# easy case: only character
axis(2, at = x_at, labels = x_labels, las = 1, cex = lab.cex,
lwd = 0, lwd.ticks = 1)
}
}
} else {
axis(2, lwd = 0, lwd.ticks = 1)
if(AXIS_AS_NUM){
axis(1, lwd = 0, lwd.ticks = 1)
} else {
if(lab.fit == "simple"){
if(any(grepl("^&", x_labels))){
# means we call expression()
# drawback => expressions can overlap
qui = grepl("^&", x_labels)
if(any(!qui)){
axis(1, at = x_at[!qui], labels = x_labels[!qui], cex.axis = lab.cex,
lwd = 0, lwd.ticks = 1)
}
for(i in which(qui)){
axis(1, at = x_at[i], labels = expr_builder(x_labels[i]), cex.axis = lab.cex,
lwd = 0, lwd.ticks = 1)
}
} else {
# easy case: only character
axis(1, at = x_at, labels = x_labels, cex.axis = lab.cex,
lwd = 0, lwd.ticks = 1)
}
} else if(lab.fit == "multi"){
axis(1, at = x_at, labels = NA, tcl = -0.25,
lwd = 0, lwd.ticks = 1)
for(my_line in unique(lab.info$line)){
qui_line = my_line == lab.info$line
x_labels_current = x_labels[qui_line]
x_at_current = x_at[qui_line]
if(any(grepl("^&", x_labels_current))){
qui = grepl("^&", x_labels_current)
if(any(!qui)){
axis(1, at = x_at_current[!qui], labels = x_labels_current[!qui],
cex.axis = lab.cex, line = my_line, lwd = 0)
}
for(i in which(qui)){
axis(1, at = x_at_current[i], labels = expr_builder(x_labels_current[i]),
cex.axis = lab.cex, line = my_line, lwd = 0)
}
} else {
# easy case: only character
axis(1, at = x_at_current, labels = x_labels_current, cex.axis = lab.cex,
line = my_line, lwd = 0, gap.axis = 0)
}
}
} else if(lab.fit == "tilted"){
axis(1, at = x_at, labels = NA, tcl = -0.25, lwd = 0, lwd.ticks = 1)
if(any(grepl("^&", x_labels))){
qui = grepl("^&", x_labels)
if(any(!qui)){
xaxis_biased(1, at = x_at[!qui], labels = x_labels[!qui], cex = lab.cex,
angle = lab.info$angle, line.min = LINE_MIN_TILTED)
}
for(i in which(qui)){
xaxis_biased(1, at = x_at[i], labels = expr_builder(x_labels[i]),
cex = lab.cex, angle = lab.info$angle, line.min = LINE_MIN_TILTED)
}
} else {
# easy case: only character
xaxis_biased(1, at = x_at, labels = x_labels, cex = lab.cex,
angle = lab.info$angle, line.min = LINE_MIN_TILTED)
}
}
}
}
}
n_id = length(unique(prms$id))
#
# pt.join ####
#
check_arg(pt.join, "logical scalar")
if(pt.join){
# We join the dots
listDefault(pt.join.par, "lwd", lwd)
listDefault(pt.join.par, "col", pt.col)
if(multiple_est){
for(i in 1:n_id){
my_join.par = pt.join.par
for(param in c("col", "lwd", "lty")){
if(!is.null(pt.join.par[[param]])){
my_join.par[[param]] = par_fit(pt.join.par[[param]], i)
}
}
my_join.par$x = prms$x[prms$id == i]
my_join.par$y = prms$y[prms$id == i]
do.call("lines", my_join.par)
}
} else {
pt.join.par$x = prms$x
pt.join.par$y = prms$y
do.call("lines", pt.join.par)
}
}
#
# The confidence intervals ####
#
x = x_value
if(!is.na(ci.lwd) && ci.lwd > 0){
# a) barre verticale
if(horiz){
segments(x0=ci_low, y0=x, x1=ci_high, y1=x, lwd = par_fit(ci.lwd, prms$id),
col = par_fit(ci.col, prms$id), lty = par_fit(ci.lty, prms$id))
} else {
segments(x0=x, y0=ci_low, x1=x, y1=ci_high, lwd = par_fit(ci.lwd, prms$id),
col = par_fit(ci.col, prms$id), lty = par_fit(ci.lty, prms$id))
}
# Formatting the bar width
if(length(ci.width) > 1){
stop("The argument 'ci.width' must be of length 1.")
}
if(is.character(ci.width)){
width_nb = tryCatch(as.numeric(gsub("%", "", ci.width)), warning = function(x) x)
if(!is.numeric(width_nb)){
stop("The value of 'ci.width' is not valid. It should be equal either to a number, either to a percentage (e.g. ci.width=\"3%\").")
}
if(grepl("%", ci.width)){
if(horiz){
total_width = diff(par("usr")[3:4])
} else {
total_width = diff(par("usr")[1:2])
}
ci.width = total_width * width_nb / 100 / 2
} else {
ci.width = width_nb / 2
}
}
# b) toppings
# Only if not a reference
qui = ci_high != ci_low
if(horiz){
# i) ci_high
segments(x0=ci_high[qui], y0=x[qui]-ci.width, x1=ci_high[qui], y1=x[qui]+ci.width,
lwd = par_fit(ci.lwd, prms$id[qui]), col = par_fit(ci.col, prms$id[qui]),
lty = par_fit(ci.lty, prms$id[qui]))
# ii) ci_low
segments(x0=ci_low[qui], y0=x[qui]-ci.width, x1=ci_low[qui], y1=x[qui]+ci.width,
lwd = par_fit(ci.lwd, prms$id[qui]), col = par_fit(ci.col, prms$id[qui]),
lty = par_fit(ci.lty, prms$id[qui]))
} else {
# i) ci_high
segments(x0=x[qui]-ci.width, y0=ci_high[qui], x1=x[qui]+ci.width, y1=ci_high[qui],
lwd = par_fit(ci.lwd, prms$id[qui]), col = par_fit(ci.col, prms$id[qui]),
lty = par_fit(ci.lty, prms$id[qui]))
# ii) ci_low
segments(x0=x[qui]-ci.width, y0=ci_low[qui], x1=x[qui]+ci.width, y1=ci_low[qui],
lwd = par_fit(ci.lwd, prms$id[qui]), col = par_fit(ci.col, prms$id[qui]),
lty = par_fit(ci.lty, prms$id[qui]))
}
}
#
# ci.join ####
#
if(ci.join){
# We join the extremities of the conf. int.
listDefault(ci.join.par, "lwd", lwd)
listDefault(ci.join.par, "col", col)
listDefault(ci.join.par, "lty", 2)
if(multiple_est){
for(i in 1:n_id){
my_join.par = ci.join.par
for(param in c("col", "lwd", "lty")){
if(!is.null(ci.join.par[[param]])){
my_join.par[[param]] = par_fit(ci.join.par[[param]], i)
}
}
if(horiz){
my_join.par$y = prms$y[prms$id == i]
my_join.par$x = prms$ci_high[prms$id == i]
do.call("lines", my_join.par)
my_join.par$x = prms$ci_low[prms$id == i]
do.call("lines", my_join.par)
} else {
my_join.par$x = prms$x[prms$id == i]
my_join.par$y = prms$ci_high[prms$id == i]
do.call("lines", my_join.par)
my_join.par$y = prms$ci_low[prms$id == i]
do.call("lines", my_join.par)
}
}
} else {
if(horiz){
ci.join.par$y = prms$y
ci.join.par$x = prms$ci_high
do.call("lines", ci.join.par)
ci.join.par$x = prms$ci_low
do.call("lines", ci.join.par)
} else {
ci.join.par$x = prms$x
ci.join.par$y = prms$ci_high
do.call("lines", ci.join.par)
ci.join.par$y = prms$ci_low
do.call("lines", ci.join.par)
}
}
}
#
# ci.fill ####
#
if(ci.fill){
# We join the extremities of the conf. int.
listDefault(ci.fill.par, "col", rgb(0.5, 0.5, 0.5))
listDefault(ci.fill.par, "alpha", 0.5)
listDefault(ci.fill.par, "border", NA)
if(multiple_est){
# We create the appropriate colors
my_cols = par_fit(ci.fill.par$col, 1:n_id)
my_alpha = par_fit(ci.fill.par$alpha, 1:n_id)
ci.fill.par$alpha = NULL
for(i in 1:n_id){
current_col = as.vector(col2rgb(my_cols[i], alpha = TRUE)) / 255
if(current_col[4] == 1){
my_cols[i] = rgb(current_col[1], current_col[2], current_col[3], my_alpha[i])
}
}
ci.fill.par$col = my_cols
for(i in 1:n_id){
my_join.par = ci.fill.par
for(param in c("col", "lwd", "lty")){
if(!is.null(ci.fill.par[[param]])){
my_join.par[[param]] = par_fit(ci.fill.par[[param]], i)
}
}
if(horiz){
my_join.par$y = c(prms$y[prms$id == i], rev(prms$y[prms$id == i]))
my_join.par$x = c(prms$ci_high[prms$id == i], rev(prms$ci_low[prms$id == i]))
do.call("polygon", my_join.par)
} else {
my_join.par$x = c(prms$x[prms$id == i], rev(prms$x[prms$id == i]))
my_join.par$y = c(prms$ci_high[prms$id == i], rev(prms$ci_low[prms$id == i]))
do.call("polygon", my_join.par)
}
}
} else {
# colors => adding alpha if needed
my_alpha = ci.fill.par$alpha[1]
ci.fill.par$alpha = NULL
my_col = as.vector(col2rgb(ci.fill.par$col[1], alpha = TRUE)) / 255
if(my_col[4] == 1){
ci.fill.par$col = rgb(my_col[1], my_col[2], my_col[3], my_alpha)
}
if(horiz){
ci.fill.par$y = c(prms$y, rev(prms$y))
ci.fill.par$x = c(prms$ci_high, rev(prms$ci_low))
do.call("polygon", ci.fill.par)
} else {
ci.fill.par$x = c(prms$x, rev(prms$x))
ci.fill.par$y = c(prms$ci_high, rev(prms$ci_low))
do.call("polygon", ci.fill.par)
}
}
}
#
# Point estimates ####
#
if(!add){
# now the points or lines
if(plot_prms$type != "n"){
point.par = plot_prms[c("x", "y", "type", "cex", "col", "pch", "lty", "lwd")]
point.par$pch = par_fit(pt.pch, prms$id)
point.par$cex = par_fit(pt.cex, prms$id)
point.par$col = par_fit(pt.col, prms$id)
point.par$lwd = par_fit(pt.lwd, prms$id)
if(!is.null(pt.bg)) point.par$bg = par_fit(pt.bg, prms$id)
point.par = point.par[lengths(point.par) > 0]
do.call("points", point.par)
}
} else {
plot_prms$pch = par_fit(pt.pch, prms$id)
plot_prms$cex = par_fit(pt.cex, prms$id)
plot_prms$col = par_fit(pt.col, prms$id)
plot_prms$lwd = par_fit(pt.lwd, prms$id)
if(!is.null(pt.bg)) point.par$bg = par_fit(pt.bg, prms$id)
do.call("points", plot_prms)
}
#
# Group ####
#
if(IS_GROUP){
if(!is.list(group)) stop("Argument 'group' must be a list.")
axis_prms_fit = c("lwd", "tcl", "line", "tick", "lty", "col", "lwd.ticks", "col.ticks")
# The line and text adjustments
line.adj = group.par$line.adj
if(is.null(line.adj)){
line.adj = 0
}
group.par$line.adj = NULL
text.adj = group.par$text.adj
if(is.null(text.adj)){
text.adj = 0
}
group.par$text.adj = NULL
# axis
listDefault(group.par, "lwd", 2)
listDefault(group.par, "tcl", 0.75)
group.par$line = group.baseline + line.adj
axis_par = group.par[!grepl("^text", names(group.par))]
axis_par$side = 1
axis_par$labels = NA
# label
listDefault(group.par, "text.line", group.par$line - 0.5 + text.adj)
listDefault(group.par, "text.cex", 1)
listDefault(group.par, "text.font", 1)
listDefault(group.par, "text.col", 1)
for(i in seq_along(group)){
group_name = names(group)[i]
my_group = group[[i]]
# formatting properly
if(is.numeric(my_group)){
g_start = min(my_group)
g_end = max(my_group)
if(any(my_group %% 1 != 0)){
stop("The elements of the argument 'group' must be integers (e.g. group=list(group_name=1:2)). Currently they are numeric but not integers. Alternatively, you could use coefficient names (see details).")
}
if(g_start < 1){
warning("The elements of the argument 'group' must be integers ranging from 1 to the number of coefficients. The value of ", g_start, " has been replaced by 1.")
}
if(g_end > nrow(prms)){
warning("The elements of the argument 'group' must be integers ranging from 1 to the number of coefficients (here ", g_end, "). The value of ", g_end, " has been replaced by ", g_end, ".")
}
} else {
if(!is.character(my_group)){
stop("The elements of the argument 'group' must be either: i) a character string of length 1 or 2, or ii) integers. Currently it is not character nor numeric.\nExample of valid use: group=list(group_name=\"pattern\"), group=list(group_name=c(\"var_start\", \"var_end\")), group=list(group_name=1:2))")
}
if(!length(my_group) %in% 1:2){
stop("The elements of the argument 'group' must be either: i) a character string of length 1 or 2, or ii) integers. Currently it is a character of length ", length(my_group), ".\nExample of valid use: group=list(group_name=\"pattern\"), group=list(group_name=c(\"var_start\", \"var_end\")), group=list(group_name=1:2))")
}
if(length(my_group) == 1){
# This is a pattern
if(grepl("^%", my_group)){
qui = grepl(gsub("^%", "", my_group), x_labels_raw)
} else {
qui = grepl(my_group, x_labels)
}
if(!any(qui)){
warning("In argument 'group', the pattern: \"", my_group, "\", did not match any coefficient name.")
next
}
} else {
# This is a pattern
check = c(FALSE, FALSE)
qui = rep(FALSE, length(x_labels))
for(i in 1:2){
val = my_group[i]
if(grepl("^%", val)){
val = gsub("^%", "", val)
check = val %in% x_labels_raw
qui = qui | x_labels_raw %in% val
} else {
check = val %in% x_labels
qui = qui | x_labels %in% val
}
}
check = my_group %in% x_labels
if(!all(check)){
warning("In argument 'group', the value", enumerate_items(my_group[check], "s.quote"), ", did not match any coefficient name.")
next
}
}
qui_nb = which(qui)
g_start = min(qui_nb)
g_end = max(qui_nb)
}
my_axis = axis_par
my_axis$at = x_at[c(g_start, g_end)]
for(p in intersect(names(my_axis), axis_prms_fit)){
my_axis[[p]] = par_fit(my_axis[[p]], i)
}
side = 1 + horiz
my_axis$side = side
info = do.call("axis", my_axis)
if(!is.null(group_name)){
if(grepl("^&", group_name)){
group_name = eval(str2lang(gsub("^&", "", group_name)))
}
axis(side, mean(info), labels = group_name, tick = FALSE,
line = par_fit(group.par$text.line, i),
cex.axis = par_fit(group.par$text.cex, i),
font = par_fit(group.par$text.font, i),
col.axis = par_fit(group.par$text.col, i))
}
}
}
res = list(prms = prms, is_iplot = is_iplot, at = x_at, labels = x_labels)
return(invisible(res))
}
coefplot_prms = function(all_models, vcov = NULL, se, ci_low, ci_high, x, x.shift = 0,
dict, i.select = NULL,
keep, drop, order, ci_level = 0.95, df.t = NULL, ref = "auto",
is_iplot = TRUE, sep, as.multiple = FALSE, is_root = TRUE){
# get the default for:
# dict, ci.level, ref
varlist = list(ci = paste0(ci_level * 100, "%"))
dots_drop = c()
suggest_ref_line = FALSE
multiple_est = FALSE
NO_NAMES = FALSE
set_up(1 + is_iplot + !is_root)
AXIS_AS_NUM = FALSE
if(is_root){
# This is a list of estimations
#
# Multiple estimations ####
#
multiple_est = TRUE
nb_est = length(all_models)
rerun = FALSE
first = TRUE
while(first || rerun){
first = FALSE
res = varlist = list()
all_inter = c()
all_inter_root = c()
dots_drop = c()
num_axis = TRUE
suggest_ref_line = TRUE
for(i in 1:nb_est){
# cat("Eval =", i, "\n")
prms = try(coefplot_prms(all_models[[i]], is_root = FALSE, vcov = vcov,
se = se, ci_low = ci_low, ci_high = ci_high, x = x,
x.shift = x.shift, dict = dict, i.select = i.select,
keep = keep, drop = drop, order = order,
ci_level = ci_level, df.t = df.t, ref = ref,
is_iplot = is_iplot, sep = sep, as.multiple = as.multiple),
silent = TRUE)
if(inherits(prms, "try-error")){
if(grepl("^[^\n]+(coefplot|iplot)", prms)){
prms = stringmagic::string_clean(prms, "^[^\n]+\n")
}
# we say if the argument is invalid
dots = get("dots", parent.frame())
nm_pblm = setdiff(names(dots), c("object", ""))
msg = ""
if(length(nm_pblm) > 0){
msg = sma("NOTA: the argument{$s, enum.bq, are ? nm_pblm} not {$(a ;)}valid argument{$s}.")
}
stop_up("The {nth ? i} model raises and error:\n{as.character(prms)}", msg)
}
if(nb_est == 1){
return(prms)
}
# Some meta variables
varlist$ci = unique(prms$varlist$ci)
varlist$depvar = unique(prms$varlist$depvar)
varlist$var = unique(prms$varlist$var)
varlist$fe = unique(prms$varlist$fe)
varlist$i = unique(prms$varlist$i)
dots_drop = unique(c(dots_drop, prms$dots_drop))
suggest_ref_line = suggest_ref_line && prms$suggest_ref_line
num_axis = num_axis && prms$num_axis
df_prms = prms$prms
df_prms$est_nb = i
res[[i]] = df_prms
}
}
AXIS_AS_NUM = num_axis
all_estimates = do.call("rbind", res)
# We respect the order provided by the user
my_names_order = unique(all_estimates$estimate_names)
my_order = 1:length(my_names_order)
names(my_order) = my_names_order
all_estimates$id_order = my_order[as.character(all_estimates$estimate_names)]
all_estimates = all_estimates[base::order(all_estimates$id_order, all_estimates$est_nb), ]
# we rescale -- beware of multiple est whose x is numeric!!!
if(nb_est > 1){
if(missing(sep)){
all_sep = c(0.2, 0.2, 0.18, 0.16)
if(length(all_sep) < nb_est - 1){
sep = 1 / (nb_est-1) * 0.7
} else {
sep = all_sep[nb_est - 1]
}
} else {
n_sep = length(sep)
if(n_sep > 1){
if(n_sep < nb_est - 1){
sep = sep[n_sep]
} else {
sep = sep[nb_est - 1]
}
}
}
if(AXIS_AS_NUM){
all_estimates$x_new = all_estimates$x + ((all_estimates$est_nb - 1) / (nb_est - 1) - 0.5) * ((nb_est - 1) * sep)
} else {
all_estimates$x_new = all_estimates$id_order + ((all_estimates$est_nb - 1) / (nb_est - 1) - 0.5) * ((nb_est - 1) * sep)
}
} else {
sep = 0
all_estimates$x_new = all_estimates$id_order
}
# The coefficients
estimate = all_estimates$y
names(estimate) = all_estimates$estimate_names
ci_low = all_estimates$ci_low
ci_high = all_estimates$ci_high
x = all_estimates$x_new
prms = data.frame(estimate = estimate, ci_low = ci_low, ci_high = ci_high,
x = x, id = all_estimates$est_nb,
estimate_names = all_estimates$estimate_names,
estimate_names_raw = all_estimates$estimate_names_raw)
if(!is.null(all_estimates$is_ref)){
prms$is_ref = all_estimates$is_ref
}
} else {
#
# Single estimation ####
#
# this is a single estimation
object = all_models
if(inherits(object, "fixest")){
mat = coeftable(object, vcov = vcov)
# special case: sunab
if(is_iplot && isTRUE(object$is_sunab)){
info_sunab = object$model_matrix_info$sunab
info_period = attr(info_sunab$agg_period, "model_matrix_info")
if(mean(info_period$coef_names_full %in% rownames(mat)) >= 0.5){
# we check that the period agg coefficients are in the coef matrix
# note that some coefs can be removed due to collin
# => we're OK
} else {
mat = coeftable(summary(object, agg = "period"), vcov = vcov)
}
object$model_matrix_info = append(list(info_period), object$model_matrix_info)
}
se = mat[, 2]
estimate = mat[, 1]
names(estimate) = rownames(mat)
if("fml" %in% names(object)){
depvar = gsub(" ", "", as.character(object$fml)[[2]])
if(depvar %in% names(dict)) depvar = dict[depvar]
varlist$depvar = depvar
}
} else if(is.list(object) && !is.null(oldClass(object))){
sum_exists = FALSE
for(c_name in class(object)){
if(exists(paste0("summary.", c_name), mode = "function")){
sum_exists = TRUE
break
}
}
# the dep var
fml = try(formula(object), silent = TRUE)
if(!inherits(fml, "try-error") && length(fml) == 3){
depvar = gsub(" ", "", as.character(fml)[[2]])
if(depvar %in% names(dict)) depvar = dict[depvar]
varlist$depvar = depvar
}
if(!sum_exists){
# stop("There is no summary method for objects of class ", c_name, ". 'coefplot' applies summary to the object to extract the coeftable. Maybe add directly the coeftable in object instead?")
mat = coeftable(object)
} else {
obj_sum = summary(object)
mat = coeftable(obj_sum)
}
m_names = tolower(colnames(mat))
if(ncol(mat) == 4 || (grepl("estimate", m_names[1]) &&
grepl("std\\.? error", m_names[1]))){
se = mat[, 2]
estimate = mat[, 1]
names(estimate) = rownames(mat)
}
} else if(is.matrix(object)){
# object is a matrix containing the coefs and SEs
m_names = tolower(colnames(object))
if(ncol(object) == 4 || (grepl("estimate", m_names[1]) &&
grepl("std\\.? error", m_names[1]))){
se = object[, 2]
estimate = object[, 1]
names(estimate) = rownames(object)
} else {
stop_up("Argument 'object' is a matrix but it should contain 4 columns (the two first ones should be reporting the estimate and the standard-error). Either provide an appropriate matrix or give directly the vector of estimated coefficients in arg. estimate.")
}
} else if(length(object[1]) > 1 || !is.null(dim(object)) || !is.numeric(object)){
stop_up("Argument 'object' must be either: i) a `fixest` estimation, ii) a matrix of coefficients table, or iii) a numeric vector of the point estimates. Currently it is neither of the three.")
} else {
# it's a numeric vector
estimate = object
}
if(!is.list(object)){
object = list()
}
n = length(estimate)
if(missing(se)){
if(missing(ci_low) || missing(ci_high)){
stop_up("When passing a vector of coefficients, the values for ",
"`se` or `ci_low` and `ci_high` must be explicitly provided.\n",
"Problem: `se`, `ci_low` and `ci_high` are missing.")
}
varlist$ci = NULL
} else {
if(!missing(ci_low) || !missing(ci_high)){
warning("Since 'se' is provided, arguments 'ci_low' or 'ci_high' are ignored.")
}
# We compute the CI
check_arg(df.t, "NULL | numeric scalar GE{1}")
if(inherits(object, "fixest")){
if(is.null(df.t)){
df.t = degrees_freedom(object, "t")
if(is.null(df.t)){
# may for user defined functions using fixest
df.t = degrees_freedom(object, "resid")
if(is.null(df.t)){
# let's be safe
df.t = Inf
}
}
}
nb = fixest_CI_factor(object, ci_level, df.t = df.t)[2]
} else {
# non fixest objects
if(is.null(df.t)){
df.t = tryCatch(nobs(object) - length(coef(object)),
error = function(e) NULL)
if(is.null(df.t)){
if(was_not_recently_used("coefplot_df")){
mema("The degrees of freedom for the t distribution could",
" not be deduced. Using a Normal distribution instead.\n",
"Note that you can provide the argument `df.t` directly.")
}
df.t = Inf
}
}
nb = abs( qt( (1-ci_level) / 2, df.t) )
}
ci_high = estimate + nb*se
ci_low = estimate - nb*se
}
#
# iplot ####
#
ref_id = NA
xlab_suggest = NULL
if(is_iplot){
if(is.null(names(estimate))){
stop_up("'iplot' must be used only with fixest objects containing variables created with i(). Currently it does not seem to be the case.")
}
all_vars = names(estimate)
if(!any(grepl("::", all_vars))){
stop_up("'iplot' must be used only with fixest objects containing variables created with i(). Currently it does not seem to be the case.")
}
# Four cases:
# - factor_var::value
# - factor_var::value:xnum
# - xnum:factor_var::value
# - factor_var::value:xfact::value
# Restriction:
# it only accepts "pure" i() variables
# We can handle only case 1, 2, and 3
# case 4 is too messy
# case 4: multiple lines + legend ?
# We take the first i() in the list
# after having applied keep_apply
# avoids bug with IVs => problem if user names the variables that way
is_IV = FALSE
if(isTRUE(object$is_iv) && identical(object$iv_stage, 2)){
all_vars = gsub("^fit_", "", all_vars)
names(estimate) = all_vars
}
all_vars = keep_apply(all_vars, keep)
mm_info = object$model_matrix_info
ANY_TWO_FACTORS = FALSE
# Finding out which to display
i_running = 0
for(i in seq_along(mm_info)){
info = mm_info[[i]]
if(isFALSE(info$is_inter_fact)){
# First case: regular variable. species::setosa
if(any(info$coef_names_full %in% all_vars)){
i_running = i_running + 1
if(i.select == i_running){
# That's the one
break
}
}
}
}
if(i_running < i.select){
# not found, maybe an interaction? second round
for(i in seq_along(mm_info)){
info = mm_info[[i]]
ANY_TWO_FACTORS = ANY_TWO_FACTORS || isTRUE(info$is_inter_fact)
if(isFALSE(info$is_inter_fact)){
# Second case: interacted variable. species::setosa:x1 ou x1:species::setosa
if(!isTRUE(info$is_inter_num)){
pattern = paste0("(:", info$f_name, "::.+)|(", info$f_name, "::[^:]+:)")
vars_inter = grep(pattern, all_vars, value = TRUE)
if(length(vars_inter) > 0){
vars_inter_unik = unique(gsub(pattern, "", vars_inter))
ok = FALSE
for(v in vars_inter_unik){
inter_before = paste0(v, ":", info$coef_names_full)
inter_after = paste0(info$coef_names_full, ":", v)
if(any(c(inter_before, inter_after) %in% all_vars)){
i_running = i_running + 1
if(i.select == i_running){
# That's the one
ok = TRUE
break
}
}
}
if(ok){
if(any(inter_before %in% all_vars)){
info$coef_names_full = inter_before
} else {
info$coef_names_full = inter_after
}
break
}
}
}
}
}
}
if(i_running < i.select){
if(i_running == 0){
if(length(mm_info) == 0){
stop_up("In iplot(), no variable created with i() found (it works only with that kind of variables)")
}
if(ANY_TWO_FACTORS){
stop_up("In iplot(), no valid variable created with i() found. Note that it does not work when i() interacts two factors.")
}
stop_up("In iplot(), no valid variable created with i() found. Note that if there are interactions, they should be created with something of the form i(factor_var, num_var), otherwise the support is limited.")
}
stop_up("In iplot(), the value of 'i.select' (=", i.select, ") exceeds the number of valid i() variables found (=", i_running, "). Please give a lower number.")
}
# "keep" here works differently => new arg. i.select?
ANY_AUTO_REF = length(info$ref_id) > 0
SHOW_REF = (identical(ref, "auto") || isTRUE(ref) || identical(ref, "all")) && ANY_AUTO_REF
SHOW_REF_FIRST = SHOW_REF && !identical(ref, "all")
# Global variables for fill_coef function
names_coef = names(estimate)
names_all = info$coef_names_full
new_names = info$items
is_rm = FALSE
ID_rm = NULL
if(ANY_AUTO_REF){
if(SHOW_REF_FIRST){
ID_rm = info$ref_id[-1]
} else if(SHOW_REF == FALSE){
ID_rm = info$ref_id
}
is_rm = length(ID_rm) > 0
}
if(is_rm){
names_all = names_all[-ID_rm]
new_names = new_names[-ID_rm]
}
fill_coef = function(coef){
# we get the vector right
res = rep(0, length(names_all))
names(res) = names_all
# we need it for CI
names(coef) = names_coef
inter_names = intersect(names_all, names_coef)
res[inter_names] = coef[inter_names]
names(res) = new_names
res
}
estimate = fill_coef(estimate)
ci_high = fill_coef(ci_high)
ci_low = fill_coef(ci_low)
estimate_names = new_names
estimate_names_raw = names_all
if(isTRUE(info$is_num) && missing(x)){
AXIS_AS_NUM = TRUE
names(estimate) = NULL
x = info$items
if(is_rm) x = x[-ID_rm]
}
# ref
if(SHOW_REF){
suggest_ref_line = length(info$ref_id) == 1 && info$is_inter_num
is_ref = seq_along(estimate) == info$ref_id[1]
} else {
is_ref = rep(FALSE, length(estimate))
}
n = length(estimate)
varlist$i = dict_apply(gsub("::.*", "", names_all[1]), dict)
}
# The DF of all the parameters
prms = data.frame(estimate = estimate, ci_low = ci_low, ci_high = ci_high)
if(is_iplot){
prms$estimate_names = estimate_names
prms$estimate_names_raw = estimate_names_raw
prms$is_ref = is_ref
} else if(!is.null(names(estimate))){
prms$estimate_names = names(estimate)
prms$estimate_names_raw = names(estimate)
} else {
NO_NAMES = TRUE
prms$estimate_names = paste0("c", 1:nrow(prms))
prms$estimate_names_raw = prms$estimate_names
}
# setting the names of the estimate
if(!missing(x)){
if(length(x) != n){
stop_up("Argument 'x' must have the same length as the number of coefficients (currently {len ? x} vs {n ? n}).")
}
if(!is.numeric(x)){
names(estimate) = x
prms$estimate_names = names(estimate)
} else if(NO_NAMES) {
AXIS_AS_NUM = TRUE
}
prms$x = x
}
# We add the reference
if(!(identical(ref, "auto") || identical(ref, "all")) && length(ref) > 0 && !isFALSE(ref)){
if(AXIS_AS_NUM){
if(!is.numeric(ref) || length(ref) > 1){
check_arg(ref, "numeric scalar")
if(is.logical(ref)){
stop_up("In this context, argument 'ref' must be a numeric scalar.")
}
} else {
names(ref) = "reference"
}
} else {
if(is.null(names(ref))){
if(!is.character(ref) || length(ref) > 1){
check_arg(ref, "character scalar", .message = "Argument 'ref' must be either: a single character, either a list or a named integer vector of length 1 (The integer gives the position of the reference among the coefficients).")
} else {
refname = ref
ref = list()
ref[[refname]] = 1
}
}
ref = unlist(ref)
if(!isScalar(ref, int = TRUE)){
reason = ifelse(length(ref) == 1, " an integer", " of length 1")
stop_up("Argument 'ref' must be either: a single character, either a list or a named integer vector of length 1. The integer gives the position of the reference among the coefficients. Currently this is not ", reason, ".")
}
}
# we recreate the parameters
n = nrow(prms)
prms$is_ref = FALSE
ref_row = data.frame(estimate = 0, ci_low = 0, ci_high = 0,
estimate_names = names(ref),
estimate_names_raw = names(ref), is_ref = TRUE)
if(AXIS_AS_NUM){
ref_row$x = unname(ref)
prms = rbind(prms, ref_row)
prms = prms[base::order(prms$x), ]
x = prms$x
} else {
prms = rbind(prms, ref_row)
if(ref > n) ref = n + 1
ids = 1:n
ids[ids >= ref] = ids[ids >= ref] + 1
prms = prms[base::order(c(ids, ref)), ]
}
}
}
n = nrow(prms)
#
# order/drop/dict ####
#
if(!is.null(prms$estimate_names)){
if(!AXIS_AS_NUM){
# dict
if(missnull(dict)){
dict = c()
} else {
prms$estimate_names = dict_apply(prms$estimate_names, dict)
}
}
# dropping some coefs
all_vars = unique(prms$estimate_names)
if(!missing(keep) && length(keep) > 0){
all_vars = keep_apply(all_vars, keep)
if(length(all_vars) == 0){
msg = ""
if(is_iplot){
msg = sma(" Didn't you mean to use 'i.select'? ",
"\nAlso note that in `iplot`, the variables names are the **values** that should be in the x-axis.",
"\nFYI, valid x-axis values are: {bq, '5|etc'K, ', 'c ? valid_vars}.",
valid_vars = unique(prms$estimate_names))
}
stop_up("Argument 'keep' has removed all variables!", msg)
}
prms = prms[prms$estimate_names %in% all_vars,]
}
if(!missing(drop) && length(drop) > 0){
all_vars = drop_apply(all_vars, drop)
if(length(all_vars) == 0){
stop_up("Argument 'drop' has removed all variables!")
}
prms = prms[prms$estimate_names %in% all_vars,]
}
if(AXIS_AS_NUM && !missing(order) && length(order) > 0){
warning("Argument 'order' is ignored since the x-axis is numeric.")
}
if(!AXIS_AS_NUM){
# ordering the coefs
if(!missing(order) && length(order) > 0){
all_vars = order_apply(all_vars, order)
my_order = 1:length(all_vars)
names(my_order) = all_vars
prms$id_order = my_order[prms$estimate_names]
prms = prms[base::order(prms$id_order), ]
}
}
estimate = prms$estimate
names(estimate) = prms$estimate_names
ci_high = prms$ci_high
ci_low = prms$ci_low
if(!is.null(prms$x)) x = prms$x
n = nrow(prms)
}
#
# we create x_labels, x_value & x_at
#
# id: used for colors/lty etc
if(!multiple_est){
if(as.multiple){
prms$id = 1:nrow(prms)
} else {
prms$id = 1
}
}
if(multiple_est){
# We don't allow x.shift
my_xlim = range(x) + c(-1, 1) * ((nb_est - 1) * sep)
x_value = x
# We allow identical aliases to display properly (they can be identified with 'group')
quoi = unique(prms[, c("estimate_names", "estimate_names_raw")])
x_labels = quoi$estimate_names
x_labels_raw = quoi$estimate_names_raw
x_at = 1:length(x_labels)
} else if(!missing(x) && is.numeric(x)){
my_xlim = range(c(x + x.shift, x - x.shift))
x_value = x + x.shift
if(NO_NAMES){
x_at = NULL
x_labels = x_labels_raw = NULL
} else {
x_at = x
x_labels = prms$estimate_names
x_labels_raw = prms$estimate_names_raw
}
} else {
x_at = 1:n
x_value = 1:n + x.shift
if(NO_NAMES){
x_labels = x_labels_raw = 1:n
} else {
x_labels = prms$estimate_names
x_labels_raw = prms$estimate_names_raw
}
my_xlim = range(c(1:n + x.shift, 1:n - x.shift)) + c(-0.5, +0.5)
}
prms$x = x_value
prms$y = prms$estimate
return(list(prms = prms, num_axis = AXIS_AS_NUM, at = x_at, labels = x_labels,
x_labels_raw = x_labels_raw, varlist = varlist,
dots_drop = dots_drop, xlim = my_xlim,
suggest_ref_line = suggest_ref_line,
multiple_est = multiple_est))
}
replace_and_make_callable = function(text, varlist, text_as_expr = FALSE){
# used to make a character string callable and substitutes variables inside text with the variables
# in varlist
# ex: "Interacted with __var__" becomes "Interacted with x_beta"
# or: "&paste(\"Interacted with \", x[beta])"
if(length(text) > 1) stop("Internal problem: length of text should not be greater than 1.")
text_split = strsplit(paste0(text, " "), "__")[[1]]
if(length(text_split) < 3){
# Nothing to be done!
return(text)
} else {
# We need to replace the variables
is_var = seq_along(text_split) %% 2 == 0
my_variables = text_split[is_var]
if(length(varlist) == 0 || any(!my_variables %in% names(varlist))){
info = "No special variable is available for this estimation."
if(length(varlist) > 0){
info = paste0("In this estimation, the only special variable", enumerate_items(paste0("__", names(varlist), "__"), "s.is.start"), ". ")
}
# warning(info, enumerate_items(paste0("__", setdiff(my_variables, names(varlist)), "__"), "is"), " not valid, thus ignored.", call. = FALSE)
return("")
not_var = !my_variables %in% names(varlist)
is_var[is_var][not_var] = FALSE
my_variables = intersect(my_variables, names(varlist))
if(length(my_variables) == 0){
return(text)
}
}
my_variables_values = varlist[my_variables]
if(any(lengths(varlist[my_variables]) > 1)){
qui = which(lengths(varlist[my_variables]) > 1)[1]
n_val = lengths(varlist[my_variables])[qui]
warning("The special variable __", my_variables[qui], "__ takes ", n_val, " values, only the first is used.", call. = FALSE)
my_variables_values = sapply(my_variables_values, function(x) x[1])
} else {
my_variables_values = unlist(my_variables_values)
}
# we prepare the result (we drop the last space)
n = length(text_split)
text_new = text_split
text_new[n] = gsub(" $", "", text_new[n])
if(nchar(text_new[n]) == 0){
text_new = text_new[-n]
is_var = is_var[-n]
}
# Do the variables contain expressions?
is_expr = grepl("^&", my_variables_values)
if(any(is_expr)){
expr_drop = function(x){
if(grepl("^&", x)){
res = gsub("^&", "", x)
if(grepl("^expression\\(", res)){
res = gsub("(^expression\\()|(\\)$)", "", res)
} else if(grepl("^substitute\\(", res)){
res = deparse(eval(str2lang(res)))
}
} else {
res = x
}
res
}
my_vars = sapply(my_variables_values, expr_drop)
if(text_as_expr){
text_new[is_var][is_expr] = my_vars[is_expr]
if(all(is_expr)){
text_new = paste0("&expression(", paste(text_new, collapse = " "), ")")
} else {
my_var_no_expr = paste0('"', my_vars, '"')[!is_expr]
new_names = paste0("x___", seq_along(my_var_no_expr))
text_new[is_var][!is_expr] = new_names
text_new = paste0("&substitute(", paste(text_new, collapse = " "), ", list(", paste0(new_names, "=", my_var_no_expr, collapse = ", "), "))")
}
} else {
text_new = paste0('"', text_new, '"')
text_new[is_var][is_expr] = my_vars[is_expr]
if(all(is_expr)){
text_new = paste0("&expression(paste(", paste(text_new, collapse = ", "), "))")
} else {
my_var_no_expr = paste0('"', my_vars, '"')[!is_expr]
new_names = paste0("x___", seq_along(my_var_no_expr))
text_new[is_var][!is_expr] = new_names
text_new = paste0("&substitute(paste(", paste(text_new, collapse = ", "), "), list(", paste0(new_names, "=", my_var_no_expr, collapse = ", "), "))")
}
}
return(text_new)
} else {
# They don't contain expressions => fine, we just replace with the variables
if(text_as_expr){
my_vars = paste0('"', my_variables_values, '"')
new_names = paste0("x___", seq_along(my_vars))
text_new[is_var] = new_names
text_new = paste0("&substitute(", paste(text_new, collapse = ""), ", list(", paste0(new_names, "=", my_vars, collapse = ", "), "))")
} else {
text_new[is_var] = my_variables_values
text_new = paste(text_new, collapse = "")
}
return(text_new)
}
}
}
expr_builder = function(x){
if(grepl("^&", x)){
my_expr = gsub("^&", "", x)
if(grepl("^(expression|substitute)\\(", my_expr)){
# direct evaluation
res = eval(str2lang(my_expr), parent.frame(2))
} else {
res = eval(str2lang(paste0("expression(", my_expr, ")")))
}
} else {
res = x
}
res
}
####
#### iplot ####
####
gen_iplot = function(){
# iplot has the same arguments as coefplot
# we make all changes in coefplot
# I automatically generate the iplot function, matching all coefplot arguments
# NOTA: I do this to avoid a discrepancy btw the current dev version
# and the package being installed
file = "./R/coefplot.R"
if(!file.exists(file)) return()
env = new.env()
source(file, local = env)
if(!exists("coefplot", envir = env)) return()
coefplot = get("coefplot", env)
coefplot_args = formals(coefplot)
arg_name = names(coefplot_args)
arg_default = sapply(coefplot_args, deparse_long)
#
# iplot
#
qui_keep = !arg_name %in% c("object", "...", "i.select", "do_iplot")
iplot_args = paste0(arg_name[qui_keep], " = ", arg_default[qui_keep], collapse = ", ")
iplot_args = gsub(" = ,", ",", iplot_args)
coefplot_call = paste0(arg_name[qui_keep], " = ", arg_name[qui_keep], collapse = ", ")
iplot_fun = paste0("iplot = function(..., i.select = 1, ", iplot_args, "){\n\n",
"\tcoefplot(..., i.select = i.select, ",
coefplot_call, ", do_iplot = TRUE)\n}")
iplot_rox = "#' @describeIn coefplot Plots the coefficients generated with i()"
# Writing the functions
intro = c("# Do not edit by hand\n# => iplot calls coefplot internally\n\n\n")
s = "\n\n\n\n"
text = c(intro, s, iplot_rox, iplot_fun, s)
update_file("R/iplot.R", text)
}
####
#### Default setting ####
####
#' Sets the defaults of coefplot
#'
#' You can set the default values of most arguments of [`coefplot`] with this function.
#'
#' @inheritParams coefplot
#'
#' @param reset Logical, default is `TRUE`. If `TRUE`, then the arguments that *are not* set during the call are reset to their "factory"-default values. If `FALSE`, on the other hand, arguments that have already been modified are not changed.
#'
#' @return
#' Doesn't return anything.
#'
#' @seealso
#' [`coefplot`]
#'
#' @examples
#'
#' # coefplot has many arguments, which makes it highly flexible.
#' # If you don't like the default style of coefplot. No worries,
#' # you can set *your* default by using the function
#' # setFixest_coefplot()
#'
#' # Estimation
#' est = feols(Petal.Length ~ Petal.Width + Sepal.Length +
#' Sepal.Width | Species, iris)
#'
#' # Plot with default style
#' coefplot(est)
#'
#' # Now we permanently change some arguments
#' dict = c("Petal.Length"="Length (Petal)", "Petal.Width"="Width (Petal)",
#' "Sepal.Length"="Length (Sepal)", "Sepal.Width"="Width (Sepal)")
#'
#' setFixest_coefplot(ci.col = 2, pt.col = "darkblue", ci.lwd = 3,
#' pt.cex = 2, pt.pch = 15, ci.width = 0, dict = dict)
#'
#' # Tadaaa!
#' coefplot(est)
#'
#' # To reset to the default settings:
#' setFixest_coefplot("all", reset = TRUE)
#' coefplot(est)
#'
setFixest_coefplot = function(style, horiz = FALSE, dict = getFixest_dict(), keep,
ci.width = "1%", ci_level = 0.95, pt.pch = 20, pt.bg = NULL,
cex = 1, pt.cex = cex, col = 1:8, pt.col = col, ci.col = col,
lwd = 1, pt.lwd = lwd, ci.lwd = lwd, ci.lty = 1, grid = TRUE,
grid.par = list(lty = 3, col = "gray"), zero = TRUE,
zero.par = list(col = "black", lwd = 1), pt.join = FALSE,
pt.join.par = list(col = pt.col, lwd = lwd), ci.join = FALSE,
ci.join.par = list(lwd = lwd, col = col, lty = 2), ci.fill = FALSE,
ci.fill.par = list(col = "lightgray", alpha = 0.5),
ref.line = "auto",
ref.line.par = list(col = "black", lty = 2), lab.cex,
lab.min.cex = 0.85,
lab.max.mar = 0.25, lab.fit = "auto", xlim.add, ylim.add, sep, bg,
group = "auto", group.par = list(lwd = 2, line = 3, tcl = 0.75),
main = "Effect on __depvar__",
value.lab = "Estimate and __ci__ Conf. Int.",
ylab = NULL, xlab = NULL, sub = NULL, reset = FALSE){
fm_cp = formals(coefplot)
arg_list = names(fm_cp)
# arg_no_default = c("object", "se", "ci_low", "ci_high", "drop", "order", "ref", "add", "only.params", "as.multiple", "...", "x", "x.shift")
# m = fm_cp[!names(fm_cp) %in% arg_no_default]
# cat(gsub(" = ,", ",", paste0(names(m), " = ", sapply(m, deparse), collapse = ", ")))
iplot_default = list()
#
# Controls
#
check_arg(style, "character scalar")
check_arg(ci.width, "scalar(numeric, character) GE{0}")
check_arg(ci_level, "numeric scalar GT{0} LT{1}")
check_arg(lwd, ci.lwd, "numeric scalar GE{0}")
check_arg(grid, zero, "logical scalar")
check_set_arg("L0 list NULL{list()}", grid.par, zero.par, pt.join.par, ref.line.par)
check_arg(reset, "logical scalar")
#
# Code
#
if(missing(style)){
style = if(reset) "all" else "default"
}
if(style == "all" && reset){
opts = list(default = list(), iplot = iplot_default)
options("fixest_coefplot" = opts)
return(invisible(NULL))
}
opts = getOption("fixest_coefplot")
if(!is.list(opts)){
warning("Wrong format of getOption('fixest_coefplot'), the options of coefplot are reset.")
opts = list(default = list(), iplot = iplot_default)
}
if(reset){
if(style == "iplot"){
my_opt = iplot_default
} else {
my_opt = list()
}
} else {
my_opt = opts[[style]]
if(is.null(my_opt)){
my_opt = list()
} else if(!is.list(opts)){
warning("Wrong format of getOption('fixest_coefplot') for style '", style, "', the options of coefplot for this style are reset.")
my_opt = list()
}
}
mc = match.call()
all_args = setdiff(names(mc), c("", "reset"))
mc = mc[all_args]
# now we find which args for which we delay evaluation
arg_var_default = lapply(fm_cp[!names(fm_cp) %in% all_args], all.vars)
arg_var_default = arg_var_default[lengths(arg_var_default) > 0]
default.eval = sapply(arg_var_default, function(x) any(x %in% all_args))
if(any(default.eval)){
default.arg = names(default.eval)[default.eval]
mc[default.arg] = fm_cp[default.arg]
# we re order
my_fact = factor(names(mc), levels = arg_list)
mc = mc[order(my_fact)]
all_args = names(mc)
}
extra_values = my_opt$extra_values
if(!is.list(extra_values)) extra_values = list()
for(arg in all_args){
my_arg = mc[[arg]]
my_arg_vars = all.vars(my_arg)
if(length(my_arg_vars) == 0 || !(any(my_arg_vars %in% setdiff(arg_list, "dict")))){
if("par" %in% all.names(my_arg)){
my_opt[[arg]] = my_arg
} else {
my_opt[[arg]] = eval(my_arg)
}
} else {
my_extra_args = setdiff(my_arg_vars, arg_list)
if(length(my_extra_args) > 0){
for(v in my_extra_args) extra_values[[v]] = eval(str2lang(v), parent.frame())
}
# control that the order is proper
arg2eval = intersect(my_arg_vars, arg_list)
order_arg = which(arg_list == arg)
order_arg2eval = sapply(arg2eval, function(x) which(arg_list == x))
if(any(order_arg2eval > order_arg)){
qui = which(order_arg2eval > order_arg)[1]
stop("In argument ", arg, ": its value (", deparse(my_arg), ") is set with the argument ", arg2eval[qui], " which occurs after. If you set an argument with the value of another argument: you must use arguments appearing before.")
}
my_opt[[arg]] = my_arg
}
}
if(length(extra_values) > 0) my_opt$extra_values = extra_values
opts[[style]] = my_opt
options("fixest_coefplot" = opts)
}
#' @rdname setFixest_coefplot
getFixest_coefplot = function(){
getOption("fixest_coefplot")
}
####
#### plot method ####
####
#' plot methods for `fixest` and `fixest_multi` objects
#'
#' Plot method reporting the coefficient estimates and their confidence intervals.
#' This is a wrapper to the more complete functions [`coefplot`] and [`iplot`].
#'
#'
#' @inheritParams coefplot
#' @inheritSection etable Arguments keep, drop and order
#'
#' @param x A `fixest` estimation, for example from [`feols`].
#' @param ... Other arguments to be passed to [`coefplot`].
#'
#' @details
#' By default `plot.fixest` runs [`coefplot`] unless the estimation includes
#' the function [`sunab`], in which case it uses [`iplot`].
#'
#' The switch to `iplot` can be made with the argument `do_iplot = TRUE`.
#'
#' @return
#' It returns invisibly the data used to create the graph.
#'
#' @examples
#'
#' #
#' # Single estimation
#' #
#'
#' est = feols(Ozone ~ Temp + Solar.R, airquality)
#' plot(est)
#'
#' # focus only on the variables
#' plot(est, drop = "Cons")
#'
#' #
#' # Multiple estimations
#' #
#'
#' est_mult = feols(Ozone ~ csw(Temp, Solar.R, Wind), airquality)
#' plot(est_mult, drop = "Const")
#'
#' #
#' # DiD estimation: Sun & Abraham
#' #
#'
#' data(base_stagg)
#' # The DiD estimation
#' res_sunab = feols(y ~ x1 + sunab(year_treated, year) | id + year, base_stagg)
#' plot(res_sunab)
#'
#'
#'
#'
plot.fixest = function(x, vcov = NULL, add = FALSE, horiz = FALSE, do_iplot = NULL,
zero = TRUE, zero.par = TRUE, dict = NULL,
keep = NULL, drop = NULL, order = NULL,
ci.width = "1%", ci_level = 0.95, plot_prms = list(),
ylim = NULL, xlim = NULL,
pch = c(20, 17, 15, 21, 24, 22), col = 1:8, cex = 1, lty = 1, lwd = 1,
pt.pch = pch, pt.bg = NULL, pt.cex = cex, pt.col = col,
ci.col = col, pt.lwd = lwd, ci.lwd = lwd, ci.lty = lty,
main = "Effect on __depvar__",
value.lab = "Estimate and __ci__ Conf. Int.",
ylab = NULL, xlab = NULL, sub = NULL, ...){
check_arg(do_iplot, "NULL logical scalar")
if(is.null(do_iplot)){
do_iplot = FALSE
if(inherits(x, "fixest") && isTRUE(x$is_sunab)){
do_iplot = TRUE
} else if(inherits(x, "fixest_multi") && isTRUE(x[[1]]$is_sunab)){
do_iplot = TRUE
}
}
coefplot(objects = x, vcov = vcov, add = add, horiz = horiz, do_iplot = do_iplot,
zero = zero, zero.par = zero.par,
dict = dict, keep = keep, drop = drop, order = order, ci.width = ci.width,
ci_level = ci_level, plot_prms = plot_prms, ylim = ylim,
xlim = xlim, pch = pch, col = col, cex = cex, lty = lty, lwd = lwd,
pt.pch = pt.pch, pt.bg = pt.bg, pt.cex = pt.cex, pt.col = pt.col,
ci.col = ci.col, pt.lwd = pt.lwd, ci.lwd = ci.lwd, ci.lty = ci.lty,
main = main, value.lab = value.lab, ylab = ylab, xlab = xlab, sub = sub,
...)
}
#' @describeIn plot.fixest
plot.fixest_multi = plot.fixest
Try the fixest package in your browser
Any scripts or data that you put into this service are public.
fixest documentation built on March 18, 2026, 9:06 a.m.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Defines functions plot.fixest getFixest_coefplot setFixest_coefplot gen_iplot expr_builder replace_and_make_callable coefplot_prms coefplot
Documented in coefplot getFixest_coefplot plot.fixest setFixest_coefplot
#----------------------------------------------#
# Author: Laurent Berge
# Date creation: Mon Feb 10 09:46:01 2020
# ~: Simple function to display results
# from multiple estimations
#----------------------------------------------#
#' Plots confidence intervals and point estimates
#'
#' This function plots the results of estimations (coefficients and confidence intervals).
#' The function `iplot` restricts the output to variables created with [`i`], either
#' interactions with factors or raw factors.
#'
#' @inheritParams etable
#' @inheritSection etable Arguments keep, drop and order
#'
#' @param ... Estimation results, which can be of the following form:
#' i) an estimation object (obtained for example from [`feols`]),
#' ii) a matrix of coefficients table, iii) a numeric vector of the point
#' estimates -- the latter requiring the extra arguments `se` or `ci_low` and `ci_high`.
#' @param objects A list of `fixest` estimation objects, or `NULL` (default). If provided,
#' the objects in `...` are ignored and the only coefficients reported are the ones in the
#' argument `objects`.
#' @param vcov Versatile argument to specify the VCOV.
#' In general, it is either a character scalar equal to a VCOV type, either a formula of the form:
#' vcov_type ~ variables. The VCOV types implemented are: "iid", "hetero" (or "HC1"),
#' "cluster", "twoway", "NW" (or "newey_west"), "DK" (or "driscoll_kraay"), and "conley".
#' It also accepts object from vcov_cluster, vcov_NW, NW, vcov_DK, DK, vcov_conley and conley.
#' It also accepts covariance matrices computed externally.
#' Finally it accepts functions to compute the covariances.
#' See the vcov documentation in the vignette.
#'
#' You can pass several VCOVs (as above) if you nest them into a list.
#' If the number of VCOVs equals the number of models, eahc VCOV is mapped to the appropriate model.
#' If there is one model and several VCOVs, or if the first element of the list is equal to
#' `"each"` or `"times"`, then the estimations will be replicated and the results
#' for each estimation and each VCOV will be reported.
#' @param se The standard errors of the estimates. It may be missing.
#' @param ci_low If `se` is not provided, the lower bound of the confidence interval.
#' For each estimate.
#' @param ci_high If `se` is not provided, the upper bound of the confidence interval.
#' For each estimate.
#' @param horiz A logical scalar, default is `FALSE`. Whether to display the confidence
#' intervals horizontally instead of vertically.
#' @param x The value of the x-axis. If missing, the names of the argument `estimate`
#' are used.
#' @param x.shift Shifts the confidence intervals bars to the left or right, depending
#' on the value of `x.shift`. Default is 0.
#' @param ci.width The width of the extremities of the confidence intervals. Default is `0.1`.
#' @param ci_level Scalar between 0 and 1: the level of the CI. By default it is equal to 0.95.
#' @param add Default is `FALSE`, if the intervals are to be added to an existing
#' graph. Note that if it is the case, then the argument `x` MUST be numeric.
#' @param xlim Numeric vector of length 2 which gives the limits of the plotting region for
#' the x-axis. The default is `NULL`, which means that it is automatically defined.
#' Use the argument `xlim.add` to simply increase or decrese the default limits.
#' @param ylim Numeric vector of length 2 which gives the limits of the plotting region for
#' the y-axis. The default is `NULL`, which means that it is automatically defined.
#' Use the argument `ylim.add` to simply increase or decrese the default limits.
#' @param pch The patch of the coefficient estimates. Default is 1 (circle).
#' This is an alias to tha argument `pt.pch`.
#' @param pt.pch The patch of the coefficient estimates. Default is 1 (circle).
#' @param cex Numeric, default is 1. Expansion factor for the points
#' @param pt.cex The size of the coefficient estimates. Default is the other argument `cex`.
#' @param col The color of the points and the confidence intervals. Default is 1
#' ("black"). Note that you can set the colors separately for each of them
#' with `pt.col` and `ci.col`.
#' @param pt.col The color of the coefficient estimates. Default is equal to the argument `col`.
#' @param ci.col The color of the confidence intervals. Default is equal to the argument `col`.
#' @param lwd General line with. Default is 1.
#' @param lty The line type of the confidence intervals. Default is 1.
#' This is an alias to the argument `ci.lty`.
#' @param pt.lwd The line width of the coefficient estimates. Default is equal to
#' the other argument `lwd`.
#' @param ci.lwd The line width of the confidence intervals. Default is equal to
#' the other argument `lwd`.
#' @param ci.lty The line type of the confidence intervals. Default is 1.
#' @param grid Logical, default is `TRUE`. Whether a grid should be displayed. You
#' can set the display of the grid with the argument `grid.par`.
#' @param grid.par List. Parameters of the grid. The default values are: `lty =
#' 3` and `col = "gray"`. You can add any graphical parameter that will be passed
#' to [`graphics::abline`]. You also have two additional arguments: use `horiz =
#' FALSE` to disable the horizontal lines, and use `vert = FALSE` to disable the
#' vertical lines. Eg: `grid.par = list(vert = FALSE, col = "red", lwd = 2)`.
#' @param zero Logical scalar, default is `TRUE`. Whether the 0 should be displayed
#' in the limits of the y-axis.
#' Note that you can set how this zero line looks like with the argument `zero.par`.
#' @param zero.par A named list of graphical parameters or a logical scalar.
#' This argument is a list containing the graphical parameters used to draw the zero-line.
#' The default value is `list(col = "black", lwd = 1)` (it's the same if `TRUE`).
#' Set it to `FALSE` to turn off the special emphasis of the zero line.
#' You can add any graphical parameter that will be passed
#' to [`graphics::abline`]. Example: `zero.par = list(col = "darkblue", lwd = 3)`.
#' @param pt.join Logical, default is `FALSE`. If `TRUE`, then the coefficient estimates
#' are joined with a line.
#' @param df.t Integer scalar or `NULL` (default). The degrees of freedom (DoF) to use
#' when computing the confidence intervals with the Student t. By default it
#' tries to capture the DoF from the estimation. To use a Normal law to compute the
#' confidence interval, use `df.t = Inf`.
#' @param pt.join.par List. Parameters of the line joining the coefficients. The
#' default values are: `col = pt.col` and `lwd = lwd`. You can add any graphical
#' parameter that will be passed to [`lines`]. Eg: `pt.join.par = list(lty = 2)`.
#' @param ref Used to add points at `y = 0` (typically to visualize reference points).
#' Either: i) "auto" (default), ii) a character vector of length 1, iii) a list
#' of length 1, iv) a named integer vector of length 1, or v) a numeric vector.
#' By default, in `iplot`, if the argument `ref` has been used in the estimation,
#' these references are automatically added. If ii), ie a character scalar, then
#' that coefficient equal to zero is added as the first coefficient. If a list or
#' a named integer vector of length 1, then the integer gives the position of the
#' reference among the coefficients and the name gives the coefficient name. A non-named
#' numeric value of `ref` only works if the x-axis is also numeric (which can happen
#' in `iplot`).
#' @param ref.line Logical or numeric, default is "auto", whose behavior depends
#' on the situation. It is `TRUE` only if: i) interactions are plotted, ii) the
#' x values are numeric and iii) a reference is found. If `TRUE`, then a vertical
#' line is drawn at the level of the reference value. Otherwise, if numeric a vertical
#' line will be drawn at that specific value.
#' @param ref.line.par List. Parameters of the vertical line on the reference. The
#' default values are: `col = "black"` and `lty = 2`. You can add any graphical
#' parameter that will be passed to [`graphics::abline`].
#' Eg: `ref.line.par = list(lty = 1, lwd = 3)`.
#' @param xlim.add A numeric vector of length 1 or 2. It represents an extension
#' factor of xlim, in percentage. Eg: `xlim.add = c(0, 0.5)` extends `xlim` of 50%
#' on the right. If of length 1, positive values represent the right, and negative
#' values the left (Eg: `xlim.add = -0.5` is equivalent to `xlim.add = c(0.5, 0)`).
#' @param ylim.add A numeric vector of length 1 or 2. It represents an extension
#' factor of ylim, in percentage. Eg: `ylim.add = c(0, 0.5)` extends `ylim` of 50%
#' on the top. If of length 1, positive values represent the top, and negative values
#' the bottom (Eg: `ylim.add = -0.5` is equivalent to `ylim.add = c(0.5, 0)`).
#' @param only.params Logical, default is `FALSE`. If `TRUE` no graphic is displayed,
#' only the values of `x` and `y` used in the plot are returned.
#' @param ... Other arguments to be passed to `summary`, if `object` is an estimation,
#' and/or to the function `plot` or `lines` (if `add = TRUE`).
#' @param sep The distance between two estimates -- only when argument `object`
#' is a list of estimation results.
#' @param as.multiple Logical: default is `FALSE`. Only when `object` is a single
#' estimation result: whether each coefficient should have a different color, line
#' type, etc. By default they all get the same style.
#' @param bg Background color for the plot. By default it is white.
#' @param ci.join Logical default to `FALSE`. Whether to join the extremities of
#' the confidence intervals. If `TRUE`, then you can set the graphical parameters
#' with the argument `ci.join.par`.
#' @param ci.join.par A list of parameters to be passed to [`graphics::lines`].
#' Only used if `ci.join=TRUE`. By default it is equal to `list(lwd = lwd, col = col, lty = 2)`.
#' @param ci.fill Logical default to `FALSE`. Whether to fill the confidence intervals
#' with a color. If `TRUE`, then you can set the graphical parameters
#' with the argument `ci.fill.par`.
#' @param ci.fill.par A list of parameters to be passed to [`graphics::polygon`].
#' Only used if `ci.fill=TRUE`. By default it is equal to `list(col = "lightgray", alpha = 0.5)`.
#' Note that `alpha` is a special parameter that adds transparency to the color
#' (ranges from 0 to 1).
#' @param group A list, default is missing. Each element of the list reports the
#' coefficients to be grouped while the name of the element is the group name. Each
#' element of the list can be either: i) a character vector of length 1, ii) of
#' length 2, or ii) a numeric vector. If equal to: i) then it is interpreted as
#' a pattern: all element fitting the regular expression will be grouped (note that
#' you can use the special character "^^" to clean the beginning of the names, see
#' example), if ii) it corresponds to the first and last elements to be grouped,
#' if iii) it corresponds to the coefficients numbers to be grouped. If equal to
#' a character vector, you can use a percentage to tell the algorithm to look at
#' the coefficients before aliasing (e.g. `"%varname"`). Example of valid uses:
#' `group=list(group_name=\"pattern\")`, `group=list(group_name=c(\"var_start\", \"var_end\"))`,
#' `group=list(group_name=1:2))`. See details.
#' @param group.par A list of parameters controlling the display of the group. The
#' parameters controlling the line are: `lwd`, `tcl` (length of the tick), `line.adj`
#' (adjustment of the position, default is 0), `tick` (whether to add the ticks),
#' `lwd.ticks`, `col.ticks`. Then the parameters controlling the text: `text.adj`
#' (adjustment of the position, default is 0), `text.cex`, `text.font`, `text.col`.
#' @param pt.bg The background color of the point estimate (when the `pt.pch` is
#' in 21 to 25). Defaults to NULL.
#' @param lab.cex The size of the labels of the coefficients. Default is missing.
#' It is automatically set by an internal algorithm which can go as low as `lab.min.cex`
#' (another argument).
#' @param lab.min.cex The minimum size of the coefficients labels, as set by the
#' internal algorithm. Default is 0.85.
#' @param lab.max.mar The maximum size the left margin can take when trying to fit
#' the coefficient labels into it (only when `horiz = TRUE`). This is used in the
#' internal algorithm fitting the coefficient labels. Default is `0.25`.
#' @param lab.fit The method to fit the coefficient labels into the plotting region
#' (only when `horiz = FALSE`). Can be `"auto"` (the default), `"simple"`, `"multi"`
#' or `"tilted"`. If `"simple"`, then the classic axis is drawn. If `"multi"`, then
#' the coefficient labels are fit horizontally across several lines, such that they
#' don't collide. If `"tilted"`, then the labels are tilted. If `"auto"`, an automatic
#' choice between the three is made.
#' @param main The title of the plot. Default is `"Effect on __depvar__"`. You can
#' use the special variable `__depvar__` to set the title (useful when you set the
#' plot default with [`setFixest_coefplot`]).
#' @param value.lab The label to appear on the side of the coefficient values. If
#' `horiz = FALSE`, the label appears in the y-axis. If `horiz = TRUE`, then it
#' appears on the x-axis. The default is equal to `"Estimate and __ci__ Conf. Int."`,
#' with `__ci__` a special variable giving the value of the confidence interval.
#' @param xlab The label of the x-axis, default is `NULL`. Note that if `horiz =
#' TRUE`, it overrides the value of the argument `value.lab`.
#' @param ylab The label of the y-axis, default is `NULL`. Note that if
#' `horiz = FALSE`, it overrides the value of the argument `value.lab`.
#' @param sub A subtitle, default is `NULL`.
#' @param style A character scalar giving the style of the plot to be used. You
#' can set styles with the function [`setFixest_coefplot`], setting all the default
#' values of the function. If missing, then it switches to either "default" or "iplot",
#' depending on the calling function.
#' @param i.select Integer scalar, default is 1. In `iplot`, used to select which
#' variable created with `i()` to select. Only used when there are several variables
#' created with [`i`]. This is an index, just try increasing numbers to hopefully
#' obtain what you want. Note that it works much better when the variables are "pure"
#' `i()` and not interacted with other variables. For example: `i(species, x1)`
#' is good while `i(species):x1` isn't. The latter will also work but the index
#' may feel weird in case there are many `i()` variables.
#' @param do_iplot Logical, default is `FALSE`. For internal use only.
#' If `TRUE`, then `iplot` is run instead of `coefplot`.
#' @param plot_prms A named list. It may contain additionnal parameters to be passed
#' to the plot.
#'
#' @seealso
#' See [`setFixest_coefplot`] to set the default values of `coefplot`, and the estimation
#' functions: e.g. [`feols`], [`fepois`][fixest::feglm], [`feglm`], [`fenegbin`][fixest::femlm].
#'
#' @section Setting custom default values:
#' The function `coefplot` dispose of many arguments to parametrize the plots. Most
#' of these arguments can be set once an for all using the function [`setFixest_coefplot`].
#' See Example 3 below for a demonstration.
#'
#' @section iplot:
#'
#' The function `iplot` restricts `coefplot` to interactions or factors created
#' with the function [`i`]. Only *one* of the i-variables will be plotted at a time.
#' If you have several i-variables, you can navigate through them with the `i.select` argument.
#'
#' The argument `i.select` is an index that will go through all the i-variables.
#' It will work well if the variables are pure, meaning not interacted with other
#' variables. If the i-variables are interacted, the index may have an odd behavior
#' but will (in most cases) work all the same, just try some numbers up until you
#' (hopefully) obtain the graph you want.
#'
#' Note, importantly, that interactions of two factor variables are (in general)
#' disregarded since they would require a 3-D plot to be properly represented.
#'
#' @section Mathematical expressions:
#'
#' You can add [`plotmath`][grDevices::plotmath] mathematical expressions in the arguments
#' `main`, `sub`, `xlab`, or `ylab`. To do so, start the character string with an ampersand.
#' For example main = "&lambda^2".
#'
#'
#' @author
#' Laurent Berge
#'
#' @examples
#'
#' #
#' # Example 1: Stacking two sets of results on the same graph
#' #
#'
#' # Estimation on Iris data with one fixed-effect (Species)
#' # + we cluster the standard-errors
#' est = feols(Petal.Length ~ Petal.Width + Sepal.Width | Species,
#' iris, vcov = "cluster")
#'
#' # Now with "regular" standard-errors
#' est_std = summary(est, vcov = "iid")
#'
#' # You can plot the two results at once
#' coefplot(est, est_std)
#'
#' # You could also use the argument vcov
#' coefplot(est, vcov = list("cluster", "iid"))
#'
#' # Alternatively, you can use the argument x.shift
#' # to do it sequentially:
#'
#' # First graph with clustered standard-errors
#' coefplot(est, x.shift = -.2)
#'
#' # 'x.shift' was used to shift the coefficients to the left.
#'
#' # Second set of results: this time with
#' # standard-errors that are not clustered.
#' coefplot(est, vcov = "iid", x.shift = .2,
#' add = TRUE, col = 2, ci.lty = 2, pch = 15)
#'
#' legend("topright", col = 1:2, pch = 20, lwd = 1, lty = 1:2,
#' legend = c("Clustered", "IID"), title = "Standard-Errors")
#'
#'
#' #
#' # Example 2: Interactions
#' #
#'
#'
#' # Now we estimate and plot the "yearly" treatment effects
#'
#' data(base_did)
#' base_inter = base_did
#'
#' # We interact the variable 'period' with the variable 'treat'
#' est_did = feols(y ~ x1 + i(period, treat, 5) | id + period, base_inter)
#'
#' # In the estimation, the variable treat is interacted
#' # with each value of period but 5, set as a reference
#'
#' # coefplot will show all the coefficients:
#' coefplot(est_did)
#'
#' # Note that the grouping of the coefficients is due to 'group = "auto"'
#'
#' # If you want to keep only the coefficients
#' # created with i() (ie the interactions), use iplot
#' iplot(est_did)
#'
#' # We can see that the graph is different from before:
#' # - only interactions are shown,
#' # - the reference is present,
#' # => this is fully flexible
#'
#' iplot(est_did, ref.line = FALSE, pt.join = TRUE)
#'
#'
#' #
#' # What if the interacted variable is not numeric?
#'
#' # Let's create a "month" variable
#' all_months = c("aug", "sept", "oct", "nov", "dec", "jan",
#' "feb", "mar", "apr", "may", "jun", "jul")
#' base_inter$period_month = all_months[base_inter$period]
#'
#' # The new estimation
#' est = feols(y ~ x1 + i(period_month, treat, "oct") | id+period, base_inter)
#' # Since 'period_month' of type character, coefplot sorts it
#' iplot(est)
#'
#' # To respect a plotting order, use a factor
#' base_inter$month_factor = factor(base_inter$period_month, levels = all_months)
#' est = feols(y ~ x1 + i(month_factor, treat, "oct") | id + period, base_inter)
#' iplot(est)
#'
#'
#' #
#' # Example 3: Setting defaults
#' #
#'
#' # coefplot has many arguments, which makes it highly flexible.
#' # If you don't like the default style of coefplot. No worries,
#' # you can set *your* default by using the function
#' # setFixest_coefplot()
#'
#' dict = c("Petal.Length"="Length (Petal)", "Petal.Width"="Width (Petal)",
#' "Sepal.Length"="Length (Sepal)", "Sepal.Width"="Width (Sepal)")
#'
#' setFixest_coefplot(ci.col = 2, pt.col = "darkblue", ci.lwd = 3,
#' pt.cex = 2, pt.pch = 15, ci.width = 0, dict = dict)
#'
#' est = feols(Petal.Length ~ Petal.Width + Sepal.Length +
#' Sepal.Width + i(Species), iris)
#'
#' # And that's it
#' coefplot(est)
#'
#' # You can set separate default values for iplot
#' setFixest_coefplot("iplot", pt.join = TRUE, pt.join.par = list(lwd = 2, lty = 2))
#' iplot(est)
#'
#' # To reset to the default settings:
#' setFixest_coefplot("all", reset = TRUE)
#' coefplot(est)
#'
#' #
#' # Example 4: group + cleaning
#' #
#'
#' # You can use the argument group to group variables
#' # You can further use the special character "^^" to clean
#' # the beginning of the coef. name: particularly useful for factors
#'
#' est = feols(Petal.Length ~ Petal.Width + Sepal.Length +
#' Sepal.Width + Species, iris)
#'
#' # No grouping:
#' coefplot(est)
#'
#' # now we group by Sepal and Species
#' coefplot(est, group = list(Sepal = "Sepal", Species = "Species"))
#'
#' # now we group + clean the beginning of the names using the special character ^^
#' coefplot(est, group = list(Sepal = "^^Sepal.", Species = "^^Species"))
#'
#'
coefplot = function(..., objects = NULL, style = NULL, se, ci_low, ci_high, df.t = NULL,
vcov = NULL, cluster = NULL,
x, x.shift = 0, horiz = FALSE,
dict = NULL, keep, drop, order, ci.width = "1%",
ci_level = 0.95, add = FALSE, plot_prms = list(),
pch = c(20, 17, 15, 21, 24, 22), col = 1:8, cex = 1, lty = 1, lwd = 1,
ylim = NULL, xlim = NULL,
pt.pch = pch,
pt.bg = NULL, pt.cex = cex, pt.col = col, ci.col = col, pt.lwd = lwd,
ci.lwd = lwd, ci.lty = lty, grid = TRUE,
grid.par = list(lty = 3, col = "gray"),
zero = TRUE, zero.par = list(col = "black", lwd = 1), pt.join = FALSE,
pt.join.par = list(col = pt.col, lwd = lwd), ci.join = FALSE,
ci.join.par = list(lwd = lwd, col = col, lty = 2), ci.fill = FALSE,
ci.fill.par = list(col = "lightgray", alpha = 0.5), ref = "auto",
ref.line = "auto", ref.line.par = list(col = "black", lty = 2), lab.cex,
lab.min.cex = 0.85, lab.max.mar = 0.25, lab.fit = "auto", xlim.add,
ylim.add, only.params = FALSE, sep, as.multiple = FALSE,
bg, group = "auto", group.par = list(lwd = 2, line = 3, tcl = 0.75),
main = "Effect on __depvar__", value.lab = "Estimate and __ci__ Conf. Int.",
ylab = NULL, xlab = NULL, sub = NULL, i.select = NULL, do_iplot = NULL){
check_set_arg(lab.fit, "match(auto, simple, multi, tilted)")
dots = list(...)
ylab_add_ci = missing(ci_low)
#
# iplot
#
# deprecated sd argument
if("sd" %in% names(dots)){
se = dots$sd
dots$sd = NULL
}
# catching the deprecated `internal.only.i` argument
if("internal.only.i" %in% names(dots)){
is_iplot = isTRUE(dots$internal.only.i)
dots$internal.only.i = NULL
} else {
is_iplot = isTRUE(do_iplot)
}
if(is_iplot){
check_arg(i.select, "NULL integer scalar GE{1}")
if(is.null(i.select)){
i.select = 1
}
}
if(!is.null(objects)){
if(!is.list(objects) || inherits(objects, "fixest")){
dots = list(objects)
} else {
dots = objects
}
}
if(length(dots) == 0){
stop_up("You must provide at least one model to plot. Currently no model is provided.",
up = 1 * is_iplot)
}
# retro-compatibility
if("object" %in% names(dots)){
all_dots = dots$object
dots_rest = dots[names(dots) != "object"]
if(!is.null(oldClass(all_dots))){
all_dots = append(list(all_dots), dots_rest)
} else {
all_dots = append(all_dots, dots_rest)
}
} else {
all_dots = dots
}
info_models = flatten_list_of_models(all_dots, accept_non_fixest = TRUE)
all_models = info_models$all_models
#
# multiple VOCVs
#
# we only keep the arg `vcov`
if(!missnull(cluster)){
vcov = oldargs_to_vcov(NULL, cluster, vcov)
}
if(is.list(vcov) && length(vcov) > 1){
# we multiply the models
is_fixest = sapply(all_models, inherits, "fixest")
if(!all(is_fixest)){
stop_up("To use the `vcov` as a list, and hence have different VCOVs ",
"for different models, you provide only `fixest` estimations.",
"\nProblem: the {nth ? which(!is_fixest)} objects are not `fixest` models.")
}
vcov_1 = vcov[[1]]
is_rep = identical(vcov_1, "times") || identical(vcov_1, "each")
if(length(all_models) > 1 && !is_rep){
stop_up("If 'vcov' is a list, it must be of the same length as the number of models, ",
"or you should add the 'each' or 'times' keyword as the first ",
"element of the list.",
"\nProblem: there are {len ? all_models} models vs {len ? vcov} elements in `vcov`.")
}
n_models = length(all_models)
all_vcov = vcov
if(is_rep || n_models == 1){
if(is_rep){
all_vcov[[1]] = NULL
}
n_vcov = length(all_vcov)
n_total = n_models * n_vcov
if(vcov_1 == "times"){
id_mod = rep(1:n_models, n_vcov)
id_vcov = rep(1:n_vcov, each = n_models)
} else {
IS_EACH = TRUE
id_mod = rep(1:n_models, each = n_vcov)
id_vcov = rep(1:n_vcov, n_models)
}
mega_models = vector("list", n_total)
for(i in 1:n_total){
mega_models[[i]] = summary(all_models[[id_mod[i]]], vcov = all_vcov[[id_vcov[i]]])
}
all_models = mega_models
} else {
for(i in seq_along(vcov)){
all_models[[i]] = summary(all_models[[i]] , vcov = all_vcov[[i]])
}
}
vcov = NULL
}
#
# Setting the default values ####
#
opts = getOption("fixest_coefplot")
check_arg(style, "NULL character scalar")
if(is.null(style)){
if(is_iplot){
style = "iplot"
} else {
style = "default"
}
}
old_style = style
style = try(match.arg(style, choices = names(opts)), silent = TRUE)
if(inherits(style, "try-error")){
warning("Unknow style '", old_style, "', using default style instead.")
style = "default"
}
if(style != "default" && !is.null(opts[["default"]])){
# there is always inheritance from the default style
my_opt = opts[["default"]]
opts_style = opts[[style]]
if(!is.null(opts_style)){
my_opt[names(opts_style)] = opts_style
}
} else {
my_opt = opts[[style]]
}
if(!is.list(my_opt)){
warning("The default values of coefplot for style '", style, "' are ill-formed and therefore reset. Use only setFixest_coefplot for setting the default values.")
setFixest_coefplot(style = style, reset = TRUE)
} else if(length(my_opt) > 0){
arg_2_add = setdiff(names(opts[["default"]]), names(my_opt))
if(length(arg_2_add) > 0){
my_opt[arg_2_add] = opts[["default"]][arg_2_add]
# we reorder
arg_list = names(formals(coefplot))
my_fact = factor(names(my_opt), levels = arg_list)
my_opt = my_opt[base::order(my_fact)]
}
if(is_iplot){
sysOrigin = sys.parent()
mc = match.call(definition = sys.function(sysOrigin), call = sys.call(sysOrigin))
} else {
mc = match.call()
}
arg2set = setdiff(names(my_opt), names(mc))
for(arg in arg2set){
my_arg = my_opt[[arg]]
if(is.name(my_arg) || is.call(my_arg)){
if(length(my_opt$extra_values) > 0){
assign(arg, eval(my_arg, my_opt$extra_values))
} else {
assign(arg, eval(my_arg))
}
} else {
assign(arg, my_arg)
}
}
}
# Set up the dictionary
dict = setup_dict(dict, check = TRUE)
if(!is.null(dict) && any(grepl("^&", names(dict)))){
# Speficic markup to identify coefplot aliases
dict_amp = dict[grepl("^&", names(dict))]
names(dict_amp) = gsub("^&", "", names(dict_amp))
dict[names(dict_amp)] = dict_amp
}
#
# Getting the parameters => function coefplot_prms
#
info = coefplot_prms(all_models = all_models, is_root = TRUE, vcov = vcov,
se = se, ci_low = ci_low, ci_high = ci_high,
x = x, x.shift = x.shift, dict = dict, keep = keep, drop = drop,
order = order, ci_level = ci_level, df.t = df.t, ref = ref,
i.select = i.select,
is_iplot = is_iplot, sep = sep, as.multiple = as.multiple)
prms = info$prms
AXIS_AS_NUM = info$num_axis
x_at = info$at
x_labels = info$labels
x_labels_raw = info$x_labels_raw
varlist = info$varlist
my_xlim = info$xlim
suggest_ref_line = info$suggest_ref_line
multiple_est = info$multiple_est
if(only.params){
return(list(prms = prms, is_iplot = is_iplot, at = x_at, labels = x_labels))
}
check_arg(plot_prms, "named list")
ci_low = prms$ci_low
ci_high = prms$ci_high
x_value = prms$x
if(horiz){
# We just reverse the x/y
tmp = prms$x
prms$x = prms$y
prms$y = tmp
}
check_arg(xlab, "NULL character vector")
if(is.null(xlab) && is_iplot){
xlab = "__i__"
}
#
# Title ####
#
check_arg(xlim, ylim, "NULL numeric vector len(2)")
# xlab / main / ylab / sub
check_arg(value.lab, main, "character scalar")
check_arg(xlab, ylab, sub, "character scalar NULL")
if(horiz){
if(is.null(xlab)) xlab = value.lab
} else {
if(is.null(ylab)) ylab = value.lab
}
if(is.null(xlab)) xlab = ""
if(is.null(ylab)) ylab = ""
if(is.null(sub)) sub = ""
main = replace_and_make_callable(main, varlist)
ylab = replace_and_make_callable(ylab, varlist)
xlab = replace_and_make_callable(xlab, varlist)
sub = replace_and_make_callable(sub, varlist)
main = expr_builder(main)
ylab = expr_builder(ylab)
xlab = expr_builder(xlab)
sub = expr_builder(sub)
plot_prms$main = ""
plot_prms$ylab = ""
plot_prms$xlab = ""
plot_prms$sub = ""
#
# group = auto + Height ####
#
# The value of group = "auto" => renaming the labels
if(identical(group, "auto") && is_iplot == FALSE){
# we change the names of interactions
qui = grepl(":", x_labels_raw)
if(any(qui)){
# we only keep the ones with ':', or '::' or '::' and ':', no more
n_colons = lengths(strsplit(x_labels_raw, ':'))
qui = qui & n_colons <= 3
all_vars_to_group = character()
if(any(qui)){
x_inter = gsub(":.+", "", x_labels_raw[qui])
tx_inter = table(x_inter)
all_vars_to_group = names(tx_inter)[tx_inter >= 2]
}
group = list()
for(i in seq_along(all_vars_to_group)){
var_left = all_vars_to_group[i]
qui_select = substr(x_labels_raw, 1, nchar(var_left)) == var_left
# we require to be next to each other
if(!all(diff(which(qui_select)) == 1)) next
x_select = x_labels_raw[qui_select]
is_inter = TRUE
n_trim = 0
group_regex = NULL
if(all(grepl("[^:]:[^:].*::", x_select))){
# treat:period::1
first_part = strsplit(x_select[1], "::")[[1]][1]
var_right = gsub(".+:", "", first_part)
n_max = nchar(first_part) + 2
} else if(all(grepl("::.*[^:]:[^:]", x_select))){
# period::1:treat
# Note that we always put 'treat' on the left
second_part = strsplit(x_select[1], "::")[[1]][2]
var_right = var_left
var_left = gsub(".+:", "", second_part)
n_max = nchar(var_right) + 2
n_trim = nchar(var_left) + 1
group_regex = paste0("%", escape_regex(var_right), "::.+:", escape_regex(var_left))
} else if(all(grepl("::", x_select))) {
is_inter = FALSE
# This is a fixest factor
n_max = nchar(var_left) + 2
} else {
# we need to find out by ourselves...
n = min(nchar(x_select[1:2]))
x_split = strsplit(substr(x_select[1:2], 1, n), "")
start_diff = which.max(x_split[[1]] != x_split[[2]])
ok = FALSE
for(n_max in n:(nchar(var_left) + 2)){
if(all(grepl(substr(x_select[1], 1, n_max), x_select, fixed = TRUE))){
ok = TRUE
break
}
}
if(!ok) next
var_right = gsub(".+:", "", substr(x_select[1], 1, n_max))
}
if(is_inter){
v_name = dict_apply(c(var_left, var_right), dict)
group_name = replace_and_make_callable("__x__ %*% (__y__ == ldots)",
list(x = v_name[1], y = v_name[2]),
text_as_expr = TRUE)
if(is.null(group_regex)){
group[[group_name]] = escape_regex(paste0("%", var_left, ":", var_right))
} else {
group[[group_name]] = group_regex
}
} else {
v_name = dict_apply(var_left, dict)
group_name = replace_and_make_callable("__x__", list(x = v_name),
text_as_expr = TRUE)
group[[group_name]] = paste0("%^", escape_regex(var_left))
}
# We update the labels
x_labels[qui_select] = substr(x_select, n_max + 1, nchar(x_select) - n_trim)
}
}
}
IS_GROUP = !identical(group, "auto") && !add && !missing(group) && !is.null(group) && length(group) > 0 && !is.null(x_labels)
line_height = par("mai")[1] / par("mar")[1]
if(IS_GROUP){
# This means there are groups!
# We compute the height of the group based on the parameters
# passed in group.par
# => it will be used to adjust the margins
# we need: tcl + text.adj
tcl = 0.75
if(!is.null(group.par$tcl)){
tcl = max(group.par$tcl)
}
line.adj = 0
if(!is.null(group.par$line.adj)){
line.adj = max(group.par$line.adj)
}
text.adj = 0
if(!is.null(group.par$text.adj)){
text.adj = max(group.par$text.adj)
}
text.cex = 1
if(!is.null(group.par$text.cex)){
text.cex = max(group.par$text.cex)
}
group.height = (0.5 + tcl + text.cex + text.adj + line.adj)
# We perform basic checks on the group
if(!is.list(group)) stop("Argument 'group' must be a list.")
# We just take care of the special group + cleaning case
# when group starts with ^^
for(i in seq_along(group)){
my_group = group[[i]]
if(! (length(my_group) == 1 && is.character(my_group) && substr(my_group, 1, 2) == "^^") ) next
if(grepl("^%", my_group)){
qui = grepl(gsub("^%", "", my_group), x_labels_raw)
} else {
qui = grepl(my_group, x_labels)
}
if(!any(qui)){
warning("In argument 'group', the pattern: \"", my_group, "\", did not match any coefficient name.")
group[[i]] = NULL
if(length(group) == 0) IS_GROUP = FALSE
next
}
group[[i]] = range(which(qui))
x_labels = gsub(my_group, "", x_labels)
}
} else {
group.height = 0
}
#
# The limits ####
#
# xlim
if(!missnull(xlim)){
check_arg(xlim, "numeric vector len(2) no na")
my_xlim = xlim
}
if(!missnull(xlim.add)){
if((!is.numeric(xlim.add) || !length(xlim.add) %in% 1:2)){
stop("Argument 'xlim.add' must be a numeric vector of length 1 or 2. It represents an extension factor of xlim, in percentage. (Eg: xlim.add = c(0, 0.5) extends xlim of 50% on the right.) If of lentgh 1, positive values represent the right, and negative values the left (Eg: xlim.add = -0.5 is equivalent to xlim.add = c(0.5, 0)).")
}
if(length(xlim.add) == 1){
if(xlim.add > 0) {
xlim.add = c(0, xlim.add)
} else {
xlim.add = c(xlim.add, 0)
}
}
x_width = diff(my_xlim)
my_xlim = my_xlim + xlim.add * x_width
}
# ylim
my_ylim = range(c(ci_low, ci_high))
if(!missnull(ylim)){
check_arg(ylim, "numeric vector len(2) no na")
my_ylim = ylim
} else if(isTRUE(zero)){
# we add 0 to xlim
if(0 > my_ylim[2]){
my_ylim[2] = 0
} else if(0 < my_ylim[1]){
my_ylim[1] = 0
}
}
if(!missnull(ylim.add)){
if((!length(ylim.add) %in% 1:2 || !is.numeric(ylim.add))){
stop("Argument 'ylim.add' must be a numeric vector of length 1 or 2. It represents an extension factor of ylim, in percentage. (Eg: ylim.add = c(0, 0.5) extends ylim of 50% on the top.) If of lentgh 1, positive values represent the top, and negative values the bottom (Eg: ylim.add = -0.5 is equivalent to ylim.add = c(0.5, 0)).")
}
if(length(ylim.add) == 1){
if(ylim.add > 0) {
ylim.add = c(0, ylim.add)
} else {
ylim.add = c(ylim.add, 0)
}
}
y_width = diff(my_ylim)
my_ylim = my_ylim + ylim.add * y_width
}
#
# Margins and cex setting ####
#
if(!missing(lab.cex)){
lab.min.cex = lab.cex
} else {
lab.cex = 1
}
if(horiz){
# we make the labels fit into the margin
xlab.line = 3
sub.line = 4
# Algorithm:
# - if lab.cex is provided:
# => we fit the labels into the margin // we extend the margin until it fits
# - if lab.cex is NOT provided:
# => we extend the margin so that the label fit. If the margin exceeds lab.max.mar% of the plot space, we reduce the cex
# etc, until it fits. If it still does not fit => we extend the margin until it fits.
#
nlines = group.height + 2
# The last 2 means 1 line on each side of the label
if(lab.cex > lab.min.cex){
# No algorithm
lab.width_in = max(strwidth(x_labels, units = "in", cex = lab.cex))
} else {
# Algorithm
max_mar_width_in = par("pin")[1] * lab.max.mar
while(nlines * line_height + lab.width_in > max_mar_width_in && lab.cex > lab.min.cex){
lab.cex = 0.95 * lab.cex
lab.width_in = max(strwidth(x_labels, units = "in", cex = lab.cex))
}
# Final step => if we go too far, we etend the margin, even beyond max_mar_width_in
if(lab.cex < lab.min.cex){
lab.cex = lab.min.cex
lab.width_in = max(strwidth(x_labels, units = "in", cex = lab.cex))
}
}
group.baseline = 2 + lab.width_in / line_height
nlines = nlines + lab.width_in / line_height
ylab.line = nlines
if(ylab != ""){
nlines = nlines + ceiling(strheight(ylab, units = "in") / line_height)
}
total_width = nlines * line_height
if(total_width > par("mai")[2]){
new_mai = par("mai")
new_mai[2] = total_width + 0.05
op = par(mai = new_mai)
on.exit(par(op))
}
my_xlim = rev(my_xlim)
} else {
# We adjust the margin only if there are groups and they
# don't fit in the original margin
# or if there is a xlab and groups at the same time
ylab.line = 3
LINE_MIN_TILTED = 0.25
LINE_MAX_TILTED = 3
if(lab.fit == "auto"){
# We want to display ALL labels
in_to_usr = diff(my_xlim) / par("pin")[1]
w_all = strwidth(x_labels, cex = lab.cex, units = "in") * in_to_usr
em = strwidth("M", units = "in") * in_to_usr
is_collided = ((w_all[-1] + w_all[-length(w_all)]) / 2 + em) > 1
while(any(is_collided) && lab.cex > lab.min.cex){
lab.cex = 0.95 * lab.cex
w_all = strwidth(x_labels, cex = lab.cex, units = "in") * in_to_usr
is_collided = ((w_all[-1] + w_all[-length(w_all)]) / 2 + em) > 1
}
# switch to multi lines
if(any(is_collided) || lab.cex < lab.min.cex){
lab.info = xaxis_labels(at = x_at, labels = x_labels, line.max = 1,
minCex = lab.min.cex, trunc = Inf, only.params = TRUE,
xlim = my_xlim)
if(lab.info$failed){
# switch to tilted
lab.fit = "tilted"
lab.info = xaxis_biased(at = x_at, labels = x_labels,
line.min = LINE_MIN_TILTED,
line.max = LINE_MAX_TILTED,
cex = seq(lab.cex, lab.min.cex, length.out = 4),
trunc = Inf, only.params = TRUE,
ylim = my_ylim)
lab.cex = lab.info$cex
nlines = lab.info$height_line + LINE_MIN_TILTED
} else {
lab.fit = "multi"
nlines = lab.info$height_line
lab.cex = lab.info$cex
}
} else {
lab.fit = "simple"
nlines = 2 * lab.cex
}
} else if(lab.fit == "multi"){
lab.info = xaxis_labels(at = x_at, labels = x_labels, line.max = 1,
minCex = lab.min.cex, trunc = Inf,
only.params = TRUE, xlim = my_xlim)
lab.cex = lab.info$cex
nlines = lab.info$height_line
if(length(unique(lab.info$line)) == 1){
lab.fit = "simple"
nlines = 2 * lab.cex
}
} else if(lab.fit == "tilted"){
lab.info = xaxis_biased(at = x_at, labels = x_labels, line.min = LINE_MIN_TILTED,
line.max = LINE_MAX_TILTED,
cex = seq(lab.cex, lab.min.cex, length.out = 4),
trunc = Inf, only.params = TRUE, ylim = my_ylim)
lab.cex = lab.info$cex
nlines = lab.info$height_line + LINE_MIN_TILTED
} else if(lab.fit == "simple"){
nlines = 2 * lab.cex
}
if(IS_GROUP){
# tcl was set before
group.baseline = nlines + 1 - 0.75 + tcl
nlines = nlines + group.height
xlab.line = nlines
} else {
xlab.line = 3
}
if(xlab != ""){
xlab.line = xlab.line + ceiling(strheight(xlab, units = "in") / line_height) - 1
}
sub.line = xlab.line + 1
if(sub != ""){
nlines = sub.line + 1
} else if(xlab != ""){
nlines = sub.line
}
total_height = nlines * line_height
if(total_height > par("mai")[1]){
new_mai = par("mai")
new_mai[1] = total_height + 0.05
op = par(mai = new_mai)
on.exit(par(op))
}
}
#
# Plot (start) ####
#
all_plot_args = unique(c(names(par()), names(formals(plot.default))))
pblm = setdiff(names(plot_prms), all_plot_args)
if(length(pblm) > 0){
plot_prms[pblm] = NULL
}
# preparation of the do.call
plot_prms$col = col
if(horiz){
plot_prms$xlim = my_ylim
plot_prms$ylim = my_xlim
} else {
plot_prms$xlim = my_xlim
plot_prms$ylim = my_ylim
}
plot_prms$x = prms$x
plot_prms$y = prms$y
plot_prms$type = "p"
#
# Plot Call ####
#
if(!add){
plot_prms$axes = FALSE
# Nude graph
first.par = plot_prms
first.par$type = "n"
do.call("plot", first.par)
# The background
if(!missing(bg)){
dx = diff(my_xlim)
dy = diff(my_ylim)
if(horiz){
rect(xleft = my_ylim[1] - dy, ybottom = my_xlim[1] - dx, xright = my_ylim[2] + dy,
ytop = my_xlim[2] + dx, col = bg)
} else {
rect(xleft = my_xlim[1] - dx, ybottom = my_ylim[1] - dy, xright = my_xlim[2] + dx,
ytop = my_ylim[2] + dy, col = bg)
}
}
# The grid
if(grid){
listDefault(grid.par, "col", "gray")
listDefault(grid.par, "lty", 3)
listDefault(grid.par, "vert", TRUE)
listDefault(grid.par, "horiz", TRUE)
vert = grid.par$vert
g.horiz = grid.par$horiz
grid.par$vert = grid.par$horiz = NULL
if(g.horiz){
do.call("hgrid", grid.par)
}
if(vert){
do.call("vgrid", grid.par)
}
}
check_arg(zero.par, "logical scalar | named list")
if(!isFALSE(zero.par)){
if(isTRUE(zero.par)){
zero.par = list(col = "black", lwd = 1)
}
listDefault(zero.par, "lwd", 1)
listDefault(zero.par, "col", "black")
if(horiz){
zero.par$v = 0
} else {
zero.par$h = 0
}
do.call("abline", zero.par)
}
# Tha annotations
if(main != ""){
title(main = main)
}
if(xlab != ""){
title(xlab = xlab, line = xlab.line)
}
if(ylab != ""){
title(ylab = ylab, line = ylab.line)
}
if(sub != ""){
title(sub = sub, line = sub.line)
}
# Reference line
check_arg(ref.line, "logical scalar | numeric vector no na | charin(auto)")
if(identical(ref.line, "auto")){
ref.line = suggest_ref_line && length(unique(prms[prms$is_ref, "estimate_names_raw"])) == 1
}
is_ref_line = !isFALSE(ref.line)
line_direction = if(horiz) "h" else "v"
if(is.logical(ref.line)){
if(ref.line) {
ref_pblm = is.null(prms$is_ref) || !any(prms$is_ref)
if(ref_pblm && !"v" %in% names(ref.line.par)){
warning("You can use the argument 'ref.line = TRUE' only when a 'natural' reference is found, or if you provided a reference with argument 'ref'. You can still draw vertical lines by using 'v' in argument 'ref.line.par'. Example: ref.line.par=list(v = ", round(x_value[floor(length(x_value)/2)]), ", col=2).")
} else if(!ref_pblm){
where = tapply(prms[prms$is_ref, "x"], prms[prms$is_ref, "estimate_names_raw"], mean)
listDefault(ref.line.par, line_direction, where)
}
}
} else {
listDefault(ref.line.par, line_direction, ref.line)
}
if(is_ref_line){
listDefault(ref.line.par, "lty", 2)
do.call("abline", ref.line.par)
}
box(bty = par("bty"))
if(horiz){
axis(1, lwd = 0, lwd.ticks = 1)
if(AXIS_AS_NUM){
axis(2, las = 1, lwd = 0, lwd.ticks = 1)
} else {
if(any(grepl("^&", x_labels))){
# means we call expression()
# drawback => expressions can overlap
qui = grepl("^&", x_labels)
if(any(!qui)){
axis(2, at = x_at[!qui], labels = x_labels[!qui], las = 1, cex = lab.cex,
lwd = 0, lwd.ticks = 1)
}
for(i in which(qui)){
axis(2, at = x_at[i], labels = expr_builder(x_labels[i]), las = 1, cex = lab.cex,
lwd = 0, lwd.ticks = 1)
}
} else {
# easy case: only character
axis(2, at = x_at, labels = x_labels, las = 1, cex = lab.cex,
lwd = 0, lwd.ticks = 1)
}
}
} else {
axis(2, lwd = 0, lwd.ticks = 1)
if(AXIS_AS_NUM){
axis(1, lwd = 0, lwd.ticks = 1)
} else {
if(lab.fit == "simple"){
if(any(grepl("^&", x_labels))){
# means we call expression()
# drawback => expressions can overlap
qui = grepl("^&", x_labels)
if(any(!qui)){
axis(1, at = x_at[!qui], labels = x_labels[!qui], cex.axis = lab.cex,
lwd = 0, lwd.ticks = 1)
}
for(i in which(qui)){
axis(1, at = x_at[i], labels = expr_builder(x_labels[i]), cex.axis = lab.cex,
lwd = 0, lwd.ticks = 1)
}
} else {
# easy case: only character
axis(1, at = x_at, labels = x_labels, cex.axis = lab.cex,
lwd = 0, lwd.ticks = 1)
}
} else if(lab.fit == "multi"){
axis(1, at = x_at, labels = NA, tcl = -0.25,
lwd = 0, lwd.ticks = 1)
for(my_line in unique(lab.info$line)){
qui_line = my_line == lab.info$line
x_labels_current = x_labels[qui_line]
x_at_current = x_at[qui_line]
if(any(grepl("^&", x_labels_current))){
qui = grepl("^&", x_labels_current)
if(any(!qui)){
axis(1, at = x_at_current[!qui], labels = x_labels_current[!qui],
cex.axis = lab.cex, line = my_line, lwd = 0)
}
for(i in which(qui)){
axis(1, at = x_at_current[i], labels = expr_builder(x_labels_current[i]),
cex.axis = lab.cex, line = my_line, lwd = 0)
}
} else {
# easy case: only character
axis(1, at = x_at_current, labels = x_labels_current, cex.axis = lab.cex,
line = my_line, lwd = 0, gap.axis = 0)
}
}
} else if(lab.fit == "tilted"){
axis(1, at = x_at, labels = NA, tcl = -0.25, lwd = 0, lwd.ticks = 1)
if(any(grepl("^&", x_labels))){
qui = grepl("^&", x_labels)
if(any(!qui)){
xaxis_biased(1, at = x_at[!qui], labels = x_labels[!qui], cex = lab.cex,
angle = lab.info$angle, line.min = LINE_MIN_TILTED)
}
for(i in which(qui)){
xaxis_biased(1, at = x_at[i], labels = expr_builder(x_labels[i]),
cex = lab.cex, angle = lab.info$angle, line.min = LINE_MIN_TILTED)
}
} else {
# easy case: only character
xaxis_biased(1, at = x_at, labels = x_labels, cex = lab.cex,
angle = lab.info$angle, line.min = LINE_MIN_TILTED)
}
}
}
}
}
n_id = length(unique(prms$id))
#
# pt.join ####
#
check_arg(pt.join, "logical scalar")
if(pt.join){
# We join the dots
listDefault(pt.join.par, "lwd", lwd)
listDefault(pt.join.par, "col", pt.col)
if(multiple_est){
for(i in 1:n_id){
my_join.par = pt.join.par
for(param in c("col", "lwd", "lty")){
if(!is.null(pt.join.par[[param]])){
my_join.par[[param]] = par_fit(pt.join.par[[param]], i)
}
}
my_join.par$x = prms$x[prms$id == i]
my_join.par$y = prms$y[prms$id == i]
do.call("lines", my_join.par)
}
} else {
pt.join.par$x = prms$x
pt.join.par$y = prms$y
do.call("lines", pt.join.par)
}
}
#
# The confidence intervals ####
#
x = x_value
if(!is.na(ci.lwd) && ci.lwd > 0){
# a) barre verticale
if(horiz){
segments(x0=ci_low, y0=x, x1=ci_high, y1=x, lwd = par_fit(ci.lwd, prms$id),
col = par_fit(ci.col, prms$id), lty = par_fit(ci.lty, prms$id))
} else {
segments(x0=x, y0=ci_low, x1=x, y1=ci_high, lwd = par_fit(ci.lwd, prms$id),
col = par_fit(ci.col, prms$id), lty = par_fit(ci.lty, prms$id))
}
# Formatting the bar width
if(length(ci.width) > 1){
stop("The argument 'ci.width' must be of length 1.")
}
if(is.character(ci.width)){
width_nb = tryCatch(as.numeric(gsub("%", "", ci.width)), warning = function(x) x)
if(!is.numeric(width_nb)){
stop("The value of 'ci.width' is not valid. It should be equal either to a number, either to a percentage (e.g. ci.width=\"3%\").")
}
if(grepl("%", ci.width)){
if(horiz){
total_width = diff(par("usr")[3:4])
} else {
total_width = diff(par("usr")[1:2])
}
ci.width = total_width * width_nb / 100 / 2
} else {
ci.width = width_nb / 2
}
}
# b) toppings
# Only if not a reference
qui = ci_high != ci_low
if(horiz){
# i) ci_high
segments(x0=ci_high[qui], y0=x[qui]-ci.width, x1=ci_high[qui], y1=x[qui]+ci.width,
lwd = par_fit(ci.lwd, prms$id[qui]), col = par_fit(ci.col, prms$id[qui]),
lty = par_fit(ci.lty, prms$id[qui]))
# ii) ci_low
segments(x0=ci_low[qui], y0=x[qui]-ci.width, x1=ci_low[qui], y1=x[qui]+ci.width,
lwd = par_fit(ci.lwd, prms$id[qui]), col = par_fit(ci.col, prms$id[qui]),
lty = par_fit(ci.lty, prms$id[qui]))
} else {
# i) ci_high
segments(x0=x[qui]-ci.width, y0=ci_high[qui], x1=x[qui]+ci.width, y1=ci_high[qui],
lwd = par_fit(ci.lwd, prms$id[qui]), col = par_fit(ci.col, prms$id[qui]),
lty = par_fit(ci.lty, prms$id[qui]))
# ii) ci_low
segments(x0=x[qui]-ci.width, y0=ci_low[qui], x1=x[qui]+ci.width, y1=ci_low[qui],
lwd = par_fit(ci.lwd, prms$id[qui]), col = par_fit(ci.col, prms$id[qui]),
lty = par_fit(ci.lty, prms$id[qui]))
}
}
#
# ci.join ####
#
if(ci.join){
# We join the extremities of the conf. int.
listDefault(ci.join.par, "lwd", lwd)
listDefault(ci.join.par, "col", col)
listDefault(ci.join.par, "lty", 2)
if(multiple_est){
for(i in 1:n_id){
my_join.par = ci.join.par
for(param in c("col", "lwd", "lty")){
if(!is.null(ci.join.par[[param]])){
my_join.par[[param]] = par_fit(ci.join.par[[param]], i)
}
}
if(horiz){
my_join.par$y = prms$y[prms$id == i]
my_join.par$x = prms$ci_high[prms$id == i]
do.call("lines", my_join.par)
my_join.par$x = prms$ci_low[prms$id == i]
do.call("lines", my_join.par)
} else {
my_join.par$x = prms$x[prms$id == i]
my_join.par$y = prms$ci_high[prms$id == i]
do.call("lines", my_join.par)
my_join.par$y = prms$ci_low[prms$id == i]
do.call("lines", my_join.par)
}
}
} else {
if(horiz){
ci.join.par$y = prms$y
ci.join.par$x = prms$ci_high
do.call("lines", ci.join.par)
ci.join.par$x = prms$ci_low
do.call("lines", ci.join.par)
} else {
ci.join.par$x = prms$x
ci.join.par$y = prms$ci_high
do.call("lines", ci.join.par)
ci.join.par$y = prms$ci_low
do.call("lines", ci.join.par)
}
}
}
#
# ci.fill ####
#
if(ci.fill){
# We join the extremities of the conf. int.
listDefault(ci.fill.par, "col", rgb(0.5, 0.5, 0.5))
listDefault(ci.fill.par, "alpha", 0.5)
listDefault(ci.fill.par, "border", NA)
if(multiple_est){
# We create the appropriate colors
my_cols = par_fit(ci.fill.par$col, 1:n_id)
my_alpha = par_fit(ci.fill.par$alpha, 1:n_id)
ci.fill.par$alpha = NULL
for(i in 1:n_id){
current_col = as.vector(col2rgb(my_cols[i], alpha = TRUE)) / 255
if(current_col[4] == 1){
my_cols[i] = rgb(current_col[1], current_col[2], current_col[3], my_alpha[i])
}
}
ci.fill.par$col = my_cols
for(i in 1:n_id){
my_join.par = ci.fill.par
for(param in c("col", "lwd", "lty")){
if(!is.null(ci.fill.par[[param]])){
my_join.par[[param]] = par_fit(ci.fill.par[[param]], i)
}
}
if(horiz){
my_join.par$y = c(prms$y[prms$id == i], rev(prms$y[prms$id == i]))
my_join.par$x = c(prms$ci_high[prms$id == i], rev(prms$ci_low[prms$id == i]))
do.call("polygon", my_join.par)
} else {
my_join.par$x = c(prms$x[prms$id == i], rev(prms$x[prms$id == i]))
my_join.par$y = c(prms$ci_high[prms$id == i], rev(prms$ci_low[prms$id == i]))
do.call("polygon", my_join.par)
}
}
} else {
# colors => adding alpha if needed
my_alpha = ci.fill.par$alpha[1]
ci.fill.par$alpha = NULL
my_col = as.vector(col2rgb(ci.fill.par$col[1], alpha = TRUE)) / 255
if(my_col[4] == 1){
ci.fill.par$col = rgb(my_col[1], my_col[2], my_col[3], my_alpha)
}
if(horiz){
ci.fill.par$y = c(prms$y, rev(prms$y))
ci.fill.par$x = c(prms$ci_high, rev(prms$ci_low))
do.call("polygon", ci.fill.par)
} else {
ci.fill.par$x = c(prms$x, rev(prms$x))
ci.fill.par$y = c(prms$ci_high, rev(prms$ci_low))
do.call("polygon", ci.fill.par)
}
}
}
#
# Point estimates ####
#
if(!add){
# now the points or lines
if(plot_prms$type != "n"){
point.par = plot_prms[c("x", "y", "type", "cex", "col", "pch", "lty", "lwd")]
point.par$pch = par_fit(pt.pch, prms$id)
point.par$cex = par_fit(pt.cex, prms$id)
point.par$col = par_fit(pt.col, prms$id)
point.par$lwd = par_fit(pt.lwd, prms$id)
if(!is.null(pt.bg)) point.par$bg = par_fit(pt.bg, prms$id)
point.par = point.par[lengths(point.par) > 0]
do.call("points", point.par)
}
} else {
plot_prms$pch = par_fit(pt.pch, prms$id)
plot_prms$cex = par_fit(pt.cex, prms$id)
plot_prms$col = par_fit(pt.col, prms$id)
plot_prms$lwd = par_fit(pt.lwd, prms$id)
if(!is.null(pt.bg)) point.par$bg = par_fit(pt.bg, prms$id)
do.call("points", plot_prms)
}
#
# Group ####
#
if(IS_GROUP){
if(!is.list(group)) stop("Argument 'group' must be a list.")
axis_prms_fit = c("lwd", "tcl", "line", "tick", "lty", "col", "lwd.ticks", "col.ticks")
# The line and text adjustments
line.adj = group.par$line.adj
if(is.null(line.adj)){
line.adj = 0
}
group.par$line.adj = NULL
text.adj = group.par$text.adj
if(is.null(text.adj)){
text.adj = 0
}
group.par$text.adj = NULL
# axis
listDefault(group.par, "lwd", 2)
listDefault(group.par, "tcl", 0.75)
group.par$line = group.baseline + line.adj
axis_par = group.par[!grepl("^text", names(group.par))]
axis_par$side = 1
axis_par$labels = NA
# label
listDefault(group.par, "text.line", group.par$line - 0.5 + text.adj)
listDefault(group.par, "text.cex", 1)
listDefault(group.par, "text.font", 1)
listDefault(group.par, "text.col", 1)
for(i in seq_along(group)){
group_name = names(group)[i]
my_group = group[[i]]
# formatting properly
if(is.numeric(my_group)){
g_start = min(my_group)
g_end = max(my_group)
if(any(my_group %% 1 != 0)){
stop("The elements of the argument 'group' must be integers (e.g. group=list(group_name=1:2)). Currently they are numeric but not integers. Alternatively, you could use coefficient names (see details).")
}
if(g_start < 1){
warning("The elements of the argument 'group' must be integers ranging from 1 to the number of coefficients. The value of ", g_start, " has been replaced by 1.")
}
if(g_end > nrow(prms)){
warning("The elements of the argument 'group' must be integers ranging from 1 to the number of coefficients (here ", g_end, "). The value of ", g_end, " has been replaced by ", g_end, ".")
}
} else {
if(!is.character(my_group)){
stop("The elements of the argument 'group' must be either: i) a character string of length 1 or 2, or ii) integers. Currently it is not character nor numeric.\nExample of valid use: group=list(group_name=\"pattern\"), group=list(group_name=c(\"var_start\", \"var_end\")), group=list(group_name=1:2))")
}
if(!length(my_group) %in% 1:2){
stop("The elements of the argument 'group' must be either: i) a character string of length 1 or 2, or ii) integers. Currently it is a character of length ", length(my_group), ".\nExample of valid use: group=list(group_name=\"pattern\"), group=list(group_name=c(\"var_start\", \"var_end\")), group=list(group_name=1:2))")
}
if(length(my_group) == 1){
# This is a pattern
if(grepl("^%", my_group)){
qui = grepl(gsub("^%", "", my_group), x_labels_raw)
} else {
qui = grepl(my_group, x_labels)
}
if(!any(qui)){
warning("In argument 'group', the pattern: \"", my_group, "\", did not match any coefficient name.")
next
}
} else {
# This is a pattern
check = c(FALSE, FALSE)
qui = rep(FALSE, length(x_labels))
for(i in 1:2){
val = my_group[i]
if(grepl("^%", val)){
val = gsub("^%", "", val)
check = val %in% x_labels_raw
qui = qui | x_labels_raw %in% val
} else {
check = val %in% x_labels
qui = qui | x_labels %in% val
}
}
check = my_group %in% x_labels
if(!all(check)){
warning("In argument 'group', the value", enumerate_items(my_group[check], "s.quote"), ", did not match any coefficient name.")
next
}
}
qui_nb = which(qui)
g_start = min(qui_nb)
g_end = max(qui_nb)
}
my_axis = axis_par
my_axis$at = x_at[c(g_start, g_end)]
for(p in intersect(names(my_axis), axis_prms_fit)){
my_axis[[p]] = par_fit(my_axis[[p]], i)
}
side = 1 + horiz
my_axis$side = side
info = do.call("axis", my_axis)
if(!is.null(group_name)){
if(grepl("^&", group_name)){
group_name = eval(str2lang(gsub("^&", "", group_name)))
}
axis(side, mean(info), labels = group_name, tick = FALSE,
line = par_fit(group.par$text.line, i),
cex.axis = par_fit(group.par$text.cex, i),
font = par_fit(group.par$text.font, i),
col.axis = par_fit(group.par$text.col, i))
}
}
}
res = list(prms = prms, is_iplot = is_iplot, at = x_at, labels = x_labels)
return(invisible(res))
}
coefplot_prms = function(all_models, vcov = NULL, se, ci_low, ci_high, x, x.shift = 0,
dict, i.select = NULL,
keep, drop, order, ci_level = 0.95, df.t = NULL, ref = "auto",
is_iplot = TRUE, sep, as.multiple = FALSE, is_root = TRUE){
# get the default for:
# dict, ci.level, ref
varlist = list(ci = paste0(ci_level * 100, "%"))
dots_drop = c()
suggest_ref_line = FALSE
multiple_est = FALSE
NO_NAMES = FALSE
set_up(1 + is_iplot + !is_root)
AXIS_AS_NUM = FALSE
if(is_root){
# This is a list of estimations
#
# Multiple estimations ####
#
multiple_est = TRUE
nb_est = length(all_models)
rerun = FALSE
first = TRUE
while(first || rerun){
first = FALSE
res = varlist = list()
all_inter = c()
all_inter_root = c()
dots_drop = c()
num_axis = TRUE
suggest_ref_line = TRUE
for(i in 1:nb_est){
# cat("Eval =", i, "\n")
prms = try(coefplot_prms(all_models[[i]], is_root = FALSE, vcov = vcov,
se = se, ci_low = ci_low, ci_high = ci_high, x = x,
x.shift = x.shift, dict = dict, i.select = i.select,
keep = keep, drop = drop, order = order,
ci_level = ci_level, df.t = df.t, ref = ref,
is_iplot = is_iplot, sep = sep, as.multiple = as.multiple),
silent = TRUE)
if(inherits(prms, "try-error")){
if(grepl("^[^\n]+(coefplot|iplot)", prms)){
prms = stringmagic::string_clean(prms, "^[^\n]+\n")
}
# we say if the argument is invalid
dots = get("dots", parent.frame())
nm_pblm = setdiff(names(dots), c("object", ""))
msg = ""
if(length(nm_pblm) > 0){
msg = sma("NOTA: the argument{$s, enum.bq, are ? nm_pblm} not {$(a ;)}valid argument{$s}.")
}
stop_up("The {nth ? i} model raises and error:\n{as.character(prms)}", msg)
}
if(nb_est == 1){
return(prms)
}
# Some meta variables
varlist$ci = unique(prms$varlist$ci)
varlist$depvar = unique(prms$varlist$depvar)
varlist$var = unique(prms$varlist$var)
varlist$fe = unique(prms$varlist$fe)
varlist$i = unique(prms$varlist$i)
dots_drop = unique(c(dots_drop, prms$dots_drop))
suggest_ref_line = suggest_ref_line && prms$suggest_ref_line
num_axis = num_axis && prms$num_axis
df_prms = prms$prms
df_prms$est_nb = i
res[[i]] = df_prms
}
}
AXIS_AS_NUM = num_axis
all_estimates = do.call("rbind", res)
# We respect the order provided by the user
my_names_order = unique(all_estimates$estimate_names)
my_order = 1:length(my_names_order)
names(my_order) = my_names_order
all_estimates$id_order = my_order[as.character(all_estimates$estimate_names)]
all_estimates = all_estimates[base::order(all_estimates$id_order, all_estimates$est_nb), ]
# we rescale -- beware of multiple est whose x is numeric!!!
if(nb_est > 1){
if(missing(sep)){
all_sep = c(0.2, 0.2, 0.18, 0.16)
if(length(all_sep) < nb_est - 1){
sep = 1 / (nb_est-1) * 0.7
} else {
sep = all_sep[nb_est - 1]
}
} else {
n_sep = length(sep)
if(n_sep > 1){
if(n_sep < nb_est - 1){
sep = sep[n_sep]
} else {
sep = sep[nb_est - 1]
}
}
}
if(AXIS_AS_NUM){
all_estimates$x_new = all_estimates$x + ((all_estimates$est_nb - 1) / (nb_est - 1) - 0.5) * ((nb_est - 1) * sep)
} else {
all_estimates$x_new = all_estimates$id_order + ((all_estimates$est_nb - 1) / (nb_est - 1) - 0.5) * ((nb_est - 1) * sep)
}
} else {
sep = 0
all_estimates$x_new = all_estimates$id_order
}
# The coefficients
estimate = all_estimates$y
names(estimate) = all_estimates$estimate_names
ci_low = all_estimates$ci_low
ci_high = all_estimates$ci_high
x = all_estimates$x_new
prms = data.frame(estimate = estimate, ci_low = ci_low, ci_high = ci_high,
x = x, id = all_estimates$est_nb,
estimate_names = all_estimates$estimate_names,
estimate_names_raw = all_estimates$estimate_names_raw)
if(!is.null(all_estimates$is_ref)){
prms$is_ref = all_estimates$is_ref
}
} else {
#
# Single estimation ####
#
# this is a single estimation
object = all_models
if(inherits(object, "fixest")){
mat = coeftable(object, vcov = vcov)
# special case: sunab
if(is_iplot && isTRUE(object$is_sunab)){
info_sunab = object$model_matrix_info$sunab
info_period = attr(info_sunab$agg_period, "model_matrix_info")
if(mean(info_period$coef_names_full %in% rownames(mat)) >= 0.5){
# we check that the period agg coefficients are in the coef matrix
# note that some coefs can be removed due to collin
# => we're OK
} else {
mat = coeftable(summary(object, agg = "period"), vcov = vcov)
}
object$model_matrix_info = append(list(info_period), object$model_matrix_info)
}
se = mat[, 2]
estimate = mat[, 1]
names(estimate) = rownames(mat)
if("fml" %in% names(object)){
depvar = gsub(" ", "", as.character(object$fml)[[2]])
if(depvar %in% names(dict)) depvar = dict[depvar]
varlist$depvar = depvar
}
} else if(is.list(object) && !is.null(oldClass(object))){
sum_exists = FALSE
for(c_name in class(object)){
if(exists(paste0("summary.", c_name), mode = "function")){
sum_exists = TRUE
break
}
}
# the dep var
fml = try(formula(object), silent = TRUE)
if(!inherits(fml, "try-error") && length(fml) == 3){
depvar = gsub(" ", "", as.character(fml)[[2]])
if(depvar %in% names(dict)) depvar = dict[depvar]
varlist$depvar = depvar
}
if(!sum_exists){
# stop("There is no summary method for objects of class ", c_name, ". 'coefplot' applies summary to the object to extract the coeftable. Maybe add directly the coeftable in object instead?")
mat = coeftable(object)
} else {
obj_sum = summary(object)
mat = coeftable(obj_sum)
}
m_names = tolower(colnames(mat))
if(ncol(mat) == 4 || (grepl("estimate", m_names[1]) &&
grepl("std\\.? error", m_names[1]))){
se = mat[, 2]
estimate = mat[, 1]
names(estimate) = rownames(mat)
}
} else if(is.matrix(object)){
# object is a matrix containing the coefs and SEs
m_names = tolower(colnames(object))
if(ncol(object) == 4 || (grepl("estimate", m_names[1]) &&
grepl("std\\.? error", m_names[1]))){
se = object[, 2]
estimate = object[, 1]
names(estimate) = rownames(object)
} else {
stop_up("Argument 'object' is a matrix but it should contain 4 columns (the two first ones should be reporting the estimate and the standard-error). Either provide an appropriate matrix or give directly the vector of estimated coefficients in arg. estimate.")
}
} else if(length(object[1]) > 1 || !is.null(dim(object)) || !is.numeric(object)){
stop_up("Argument 'object' must be either: i) a `fixest` estimation, ii) a matrix of coefficients table, or iii) a numeric vector of the point estimates. Currently it is neither of the three.")
} else {
# it's a numeric vector
estimate = object
}
if(!is.list(object)){
object = list()
}
n = length(estimate)
if(missing(se)){
if(missing(ci_low) || missing(ci_high)){
stop_up("When passing a vector of coefficients, the values for ",
"`se` or `ci_low` and `ci_high` must be explicitly provided.\n",
"Problem: `se`, `ci_low` and `ci_high` are missing.")
}
varlist$ci = NULL
} else {
if(!missing(ci_low) || !missing(ci_high)){
warning("Since 'se' is provided, arguments 'ci_low' or 'ci_high' are ignored.")
}
# We compute the CI
check_arg(df.t, "NULL | numeric scalar GE{1}")
if(inherits(object, "fixest")){
if(is.null(df.t)){
df.t = degrees_freedom(object, "t")
if(is.null(df.t)){
# may for user defined functions using fixest
df.t = degrees_freedom(object, "resid")
if(is.null(df.t)){
# let's be safe
df.t = Inf
}
}
}
nb = fixest_CI_factor(object, ci_level, df.t = df.t)[2]
} else {
# non fixest objects
if(is.null(df.t)){
df.t = tryCatch(nobs(object) - length(coef(object)),
error = function(e) NULL)
if(is.null(df.t)){
if(was_not_recently_used("coefplot_df")){
mema("The degrees of freedom for the t distribution could",
" not be deduced. Using a Normal distribution instead.\n",
"Note that you can provide the argument `df.t` directly.")
}
df.t = Inf
}
}
nb = abs( qt( (1-ci_level) / 2, df.t) )
}
ci_high = estimate + nb*se
ci_low = estimate - nb*se
}
#
# iplot ####
#
ref_id = NA
xlab_suggest = NULL
if(is_iplot){
if(is.null(names(estimate))){
stop_up("'iplot' must be used only with fixest objects containing variables created with i(). Currently it does not seem to be the case.")
}
all_vars = names(estimate)
if(!any(grepl("::", all_vars))){
stop_up("'iplot' must be used only with fixest objects containing variables created with i(). Currently it does not seem to be the case.")
}
# Four cases:
# - factor_var::value
# - factor_var::value:xnum
# - xnum:factor_var::value
# - factor_var::value:xfact::value
# Restriction:
# it only accepts "pure" i() variables
# We can handle only case 1, 2, and 3
# case 4 is too messy
# case 4: multiple lines + legend ?
# We take the first i() in the list
# after having applied keep_apply
# avoids bug with IVs => problem if user names the variables that way
is_IV = FALSE
if(isTRUE(object$is_iv) && identical(object$iv_stage, 2)){
all_vars = gsub("^fit_", "", all_vars)
names(estimate) = all_vars
}
all_vars = keep_apply(all_vars, keep)
mm_info = object$model_matrix_info
ANY_TWO_FACTORS = FALSE
# Finding out which to display
i_running = 0
for(i in seq_along(mm_info)){
info = mm_info[[i]]
if(isFALSE(info$is_inter_fact)){
# First case: regular variable. species::setosa
if(any(info$coef_names_full %in% all_vars)){
i_running = i_running + 1
if(i.select == i_running){
# That's the one
break
}
}
}
}
if(i_running < i.select){
# not found, maybe an interaction? second round
for(i in seq_along(mm_info)){
info = mm_info[[i]]
ANY_TWO_FACTORS = ANY_TWO_FACTORS || isTRUE(info$is_inter_fact)
if(isFALSE(info$is_inter_fact)){
# Second case: interacted variable. species::setosa:x1 ou x1:species::setosa
if(!isTRUE(info$is_inter_num)){
pattern = paste0("(:", info$f_name, "::.+)|(", info$f_name, "::[^:]+:)")
vars_inter = grep(pattern, all_vars, value = TRUE)
if(length(vars_inter) > 0){
vars_inter_unik = unique(gsub(pattern, "", vars_inter))
ok = FALSE
for(v in vars_inter_unik){
inter_before = paste0(v, ":", info$coef_names_full)
inter_after = paste0(info$coef_names_full, ":", v)
if(any(c(inter_before, inter_after) %in% all_vars)){
i_running = i_running + 1
if(i.select == i_running){
# That's the one
ok = TRUE
break
}
}
}
if(ok){
if(any(inter_before %in% all_vars)){
info$coef_names_full = inter_before
} else {
info$coef_names_full = inter_after
}
break
}
}
}
}
}
}
if(i_running < i.select){
if(i_running == 0){
if(length(mm_info) == 0){
stop_up("In iplot(), no variable created with i() found (it works only with that kind of variables)")
}
if(ANY_TWO_FACTORS){
stop_up("In iplot(), no valid variable created with i() found. Note that it does not work when i() interacts two factors.")
}
stop_up("In iplot(), no valid variable created with i() found. Note that if there are interactions, they should be created with something of the form i(factor_var, num_var), otherwise the support is limited.")
}
stop_up("In iplot(), the value of 'i.select' (=", i.select, ") exceeds the number of valid i() variables found (=", i_running, "). Please give a lower number.")
}
# "keep" here works differently => new arg. i.select?
ANY_AUTO_REF = length(info$ref_id) > 0
SHOW_REF = (identical(ref, "auto") || isTRUE(ref) || identical(ref, "all")) && ANY_AUTO_REF
SHOW_REF_FIRST = SHOW_REF && !identical(ref, "all")
# Global variables for fill_coef function
names_coef = names(estimate)
names_all = info$coef_names_full
new_names = info$items
is_rm = FALSE
ID_rm = NULL
if(ANY_AUTO_REF){
if(SHOW_REF_FIRST){
ID_rm = info$ref_id[-1]
} else if(SHOW_REF == FALSE){
ID_rm = info$ref_id
}
is_rm = length(ID_rm) > 0
}
if(is_rm){
names_all = names_all[-ID_rm]
new_names = new_names[-ID_rm]
}
fill_coef = function(coef){
# we get the vector right
res = rep(0, length(names_all))
names(res) = names_all
# we need it for CI
names(coef) = names_coef
inter_names = intersect(names_all, names_coef)
res[inter_names] = coef[inter_names]
names(res) = new_names
res
}
estimate = fill_coef(estimate)
ci_high = fill_coef(ci_high)
ci_low = fill_coef(ci_low)
estimate_names = new_names
estimate_names_raw = names_all
if(isTRUE(info$is_num) && missing(x)){
AXIS_AS_NUM = TRUE
names(estimate) = NULL
x = info$items
if(is_rm) x = x[-ID_rm]
}
# ref
if(SHOW_REF){
suggest_ref_line = length(info$ref_id) == 1 && info$is_inter_num
is_ref = seq_along(estimate) == info$ref_id[1]
} else {
is_ref = rep(FALSE, length(estimate))
}
n = length(estimate)
varlist$i = dict_apply(gsub("::.*", "", names_all[1]), dict)
}
# The DF of all the parameters
prms = data.frame(estimate = estimate, ci_low = ci_low, ci_high = ci_high)
if(is_iplot){
prms$estimate_names = estimate_names
prms$estimate_names_raw = estimate_names_raw
prms$is_ref = is_ref
} else if(!is.null(names(estimate))){
prms$estimate_names = names(estimate)
prms$estimate_names_raw = names(estimate)
} else {
NO_NAMES = TRUE
prms$estimate_names = paste0("c", 1:nrow(prms))
prms$estimate_names_raw = prms$estimate_names
}
# setting the names of the estimate
if(!missing(x)){
if(length(x) != n){
stop_up("Argument 'x' must have the same length as the number of coefficients (currently {len ? x} vs {n ? n}).")
}
if(!is.numeric(x)){
names(estimate) = x
prms$estimate_names = names(estimate)
} else if(NO_NAMES) {
AXIS_AS_NUM = TRUE
}
prms$x = x
}
# We add the reference
if(!(identical(ref, "auto") || identical(ref, "all")) && length(ref) > 0 && !isFALSE(ref)){
if(AXIS_AS_NUM){
if(!is.numeric(ref) || length(ref) > 1){
check_arg(ref, "numeric scalar")
if(is.logical(ref)){
stop_up("In this context, argument 'ref' must be a numeric scalar.")
}
} else {
names(ref) = "reference"
}
} else {
if(is.null(names(ref))){
if(!is.character(ref) || length(ref) > 1){
check_arg(ref, "character scalar", .message = "Argument 'ref' must be either: a single character, either a list or a named integer vector of length 1 (The integer gives the position of the reference among the coefficients).")
} else {
refname = ref
ref = list()
ref[[refname]] = 1
}
}
ref = unlist(ref)
if(!isScalar(ref, int = TRUE)){
reason = ifelse(length(ref) == 1, " an integer", " of length 1")
stop_up("Argument 'ref' must be either: a single character, either a list or a named integer vector of length 1. The integer gives the position of the reference among the coefficients. Currently this is not ", reason, ".")
}
}
# we recreate the parameters
n = nrow(prms)
prms$is_ref = FALSE
ref_row = data.frame(estimate = 0, ci_low = 0, ci_high = 0,
estimate_names = names(ref),
estimate_names_raw = names(ref), is_ref = TRUE)
if(AXIS_AS_NUM){
ref_row$x = unname(ref)
prms = rbind(prms, ref_row)
prms = prms[base::order(prms$x), ]
x = prms$x
} else {
prms = rbind(prms, ref_row)
if(ref > n) ref = n + 1
ids = 1:n
ids[ids >= ref] = ids[ids >= ref] + 1
prms = prms[base::order(c(ids, ref)), ]
}
}
}
n = nrow(prms)
#
# order/drop/dict ####
#
if(!is.null(prms$estimate_names)){
if(!AXIS_AS_NUM){
# dict
if(missnull(dict)){
dict = c()
} else {
prms$estimate_names = dict_apply(prms$estimate_names, dict)
}
}
# dropping some coefs
all_vars = unique(prms$estimate_names)
if(!missing(keep) && length(keep) > 0){
all_vars = keep_apply(all_vars, keep)
if(length(all_vars) == 0){
msg = ""
if(is_iplot){
msg = sma(" Didn't you mean to use 'i.select'? ",
"\nAlso note that in `iplot`, the variables names are the **values** that should be in the x-axis.",
"\nFYI, valid x-axis values are: {bq, '5|etc'K, ', 'c ? valid_vars}.",
valid_vars = unique(prms$estimate_names))
}
stop_up("Argument 'keep' has removed all variables!", msg)
}
prms = prms[prms$estimate_names %in% all_vars,]
}
if(!missing(drop) && length(drop) > 0){
all_vars = drop_apply(all_vars, drop)
if(length(all_vars) == 0){
stop_up("Argument 'drop' has removed all variables!")
}
prms = prms[prms$estimate_names %in% all_vars,]
}
if(AXIS_AS_NUM && !missing(order) && length(order) > 0){
warning("Argument 'order' is ignored since the x-axis is numeric.")
}
if(!AXIS_AS_NUM){
# ordering the coefs
if(!missing(order) && length(order) > 0){
all_vars = order_apply(all_vars, order)
my_order = 1:length(all_vars)
names(my_order) = all_vars
prms$id_order = my_order[prms$estimate_names]
prms = prms[base::order(prms$id_order), ]
}
}
estimate = prms$estimate
names(estimate) = prms$estimate_names
ci_high = prms$ci_high
ci_low = prms$ci_low
if(!is.null(prms$x)) x = prms$x
n = nrow(prms)
}
#
# we create x_labels, x_value & x_at
#
# id: used for colors/lty etc
if(!multiple_est){
if(as.multiple){
prms$id = 1:nrow(prms)
} else {
prms$id = 1
}
}
if(multiple_est){
# We don't allow x.shift
my_xlim = range(x) + c(-1, 1) * ((nb_est - 1) * sep)
x_value = x
# We allow identical aliases to display properly (they can be identified with 'group')
quoi = unique(prms[, c("estimate_names", "estimate_names_raw")])
x_labels = quoi$estimate_names
x_labels_raw = quoi$estimate_names_raw
x_at = 1:length(x_labels)
} else if(!missing(x) && is.numeric(x)){
my_xlim = range(c(x + x.shift, x - x.shift))
x_value = x + x.shift
if(NO_NAMES){
x_at = NULL
x_labels = x_labels_raw = NULL
} else {
x_at = x
x_labels = prms$estimate_names
x_labels_raw = prms$estimate_names_raw
}
} else {
x_at = 1:n
x_value = 1:n + x.shift
if(NO_NAMES){
x_labels = x_labels_raw = 1:n
} else {
x_labels = prms$estimate_names
x_labels_raw = prms$estimate_names_raw
}
my_xlim = range(c(1:n + x.shift, 1:n - x.shift)) + c(-0.5, +0.5)
}
prms$x = x_value
prms$y = prms$estimate
return(list(prms = prms, num_axis = AXIS_AS_NUM, at = x_at, labels = x_labels,
x_labels_raw = x_labels_raw, varlist = varlist,
dots_drop = dots_drop, xlim = my_xlim,
suggest_ref_line = suggest_ref_line,
multiple_est = multiple_est))
}
replace_and_make_callable = function(text, varlist, text_as_expr = FALSE){
# used to make a character string callable and substitutes variables inside text with the variables
# in varlist
# ex: "Interacted with __var__" becomes "Interacted with x_beta"
# or: "&paste(\"Interacted with \", x[beta])"
if(length(text) > 1) stop("Internal problem: length of text should not be greater than 1.")
text_split = strsplit(paste0(text, " "), "__")[[1]]
if(length(text_split) < 3){
# Nothing to be done!
return(text)
} else {
# We need to replace the variables
is_var = seq_along(text_split) %% 2 == 0
my_variables = text_split[is_var]
if(length(varlist) == 0 || any(!my_variables %in% names(varlist))){
info = "No special variable is available for this estimation."
if(length(varlist) > 0){
info = paste0("In this estimation, the only special variable", enumerate_items(paste0("__", names(varlist), "__"), "s.is.start"), ". ")
}
# warning(info, enumerate_items(paste0("__", setdiff(my_variables, names(varlist)), "__"), "is"), " not valid, thus ignored.", call. = FALSE)
return("")
not_var = !my_variables %in% names(varlist)
is_var[is_var][not_var] = FALSE
my_variables = intersect(my_variables, names(varlist))
if(length(my_variables) == 0){
return(text)
}
}
my_variables_values = varlist[my_variables]
if(any(lengths(varlist[my_variables]) > 1)){
qui = which(lengths(varlist[my_variables]) > 1)[1]
n_val = lengths(varlist[my_variables])[qui]
warning("The special variable __", my_variables[qui], "__ takes ", n_val, " values, only the first is used.", call. = FALSE)
my_variables_values = sapply(my_variables_values, function(x) x[1])
} else {
my_variables_values = unlist(my_variables_values)
}
# we prepare the result (we drop the last space)
n = length(text_split)
text_new = text_split
text_new[n] = gsub(" $", "", text_new[n])
if(nchar(text_new[n]) == 0){
text_new = text_new[-n]
is_var = is_var[-n]
}
# Do the variables contain expressions?
is_expr = grepl("^&", my_variables_values)
if(any(is_expr)){
expr_drop = function(x){
if(grepl("^&", x)){
res = gsub("^&", "", x)
if(grepl("^expression\\(", res)){
res = gsub("(^expression\\()|(\\)$)", "", res)
} else if(grepl("^substitute\\(", res)){
res = deparse(eval(str2lang(res)))
}
} else {
res = x
}
res
}
my_vars = sapply(my_variables_values, expr_drop)
if(text_as_expr){
text_new[is_var][is_expr] = my_vars[is_expr]
if(all(is_expr)){
text_new = paste0("&expression(", paste(text_new, collapse = " "), ")")
} else {
my_var_no_expr = paste0('"', my_vars, '"')[!is_expr]
new_names = paste0("x___", seq_along(my_var_no_expr))
text_new[is_var][!is_expr] = new_names
text_new = paste0("&substitute(", paste(text_new, collapse = " "), ", list(", paste0(new_names, "=", my_var_no_expr, collapse = ", "), "))")
}
} else {
text_new = paste0('"', text_new, '"')
text_new[is_var][is_expr] = my_vars[is_expr]
if(all(is_expr)){
text_new = paste0("&expression(paste(", paste(text_new, collapse = ", "), "))")
} else {
my_var_no_expr = paste0('"', my_vars, '"')[!is_expr]
new_names = paste0("x___", seq_along(my_var_no_expr))
text_new[is_var][!is_expr] = new_names
text_new = paste0("&substitute(paste(", paste(text_new, collapse = ", "), "), list(", paste0(new_names, "=", my_var_no_expr, collapse = ", "), "))")
}
}
return(text_new)
} else {
# They don't contain expressions => fine, we just replace with the variables
if(text_as_expr){
my_vars = paste0('"', my_variables_values, '"')
new_names = paste0("x___", seq_along(my_vars))
text_new[is_var] = new_names
text_new = paste0("&substitute(", paste(text_new, collapse = ""), ", list(", paste0(new_names, "=", my_vars, collapse = ", "), "))")
} else {
text_new[is_var] = my_variables_values
text_new = paste(text_new, collapse = "")
}
return(text_new)
}
}
}
expr_builder = function(x){
if(grepl("^&", x)){
my_expr = gsub("^&", "", x)
if(grepl("^(expression|substitute)\\(", my_expr)){
# direct evaluation
res = eval(str2lang(my_expr), parent.frame(2))
} else {
res = eval(str2lang(paste0("expression(", my_expr, ")")))
}
} else {
res = x
}
res
}
####
#### iplot ####
####
gen_iplot = function(){
# iplot has the same arguments as coefplot
# we make all changes in coefplot
# I automatically generate the iplot function, matching all coefplot arguments
# NOTA: I do this to avoid a discrepancy btw the current dev version
# and the package being installed
file = "./R/coefplot.R"
if(!file.exists(file)) return()
env = new.env()
source(file, local = env)
if(!exists("coefplot", envir = env)) return()
coefplot = get("coefplot", env)
coefplot_args = formals(coefplot)
arg_name = names(coefplot_args)
arg_default = sapply(coefplot_args, deparse_long)
#
# iplot
#
qui_keep = !arg_name %in% c("object", "...", "i.select", "do_iplot")
iplot_args = paste0(arg_name[qui_keep], " = ", arg_default[qui_keep], collapse = ", ")
iplot_args = gsub(" = ,", ",", iplot_args)
coefplot_call = paste0(arg_name[qui_keep], " = ", arg_name[qui_keep], collapse = ", ")
iplot_fun = paste0("iplot = function(..., i.select = 1, ", iplot_args, "){\n\n",
"\tcoefplot(..., i.select = i.select, ",
coefplot_call, ", do_iplot = TRUE)\n}")
iplot_rox = "#' @describeIn coefplot Plots the coefficients generated with i()"
# Writing the functions
intro = c("# Do not edit by hand\n# => iplot calls coefplot internally\n\n\n")
s = "\n\n\n\n"
text = c(intro, s, iplot_rox, iplot_fun, s)
update_file("R/iplot.R", text)
}
####
#### Default setting ####
####
#' Sets the defaults of coefplot
#'
#' You can set the default values of most arguments of [`coefplot`] with this function.
#'
#' @inheritParams coefplot
#'
#' @param reset Logical, default is `TRUE`. If `TRUE`, then the arguments that *are not* set during the call are reset to their "factory"-default values. If `FALSE`, on the other hand, arguments that have already been modified are not changed.
#'
#' @return
#' Doesn't return anything.
#'
#' @seealso
#' [`coefplot`]
#'
#' @examples
#'
#' # coefplot has many arguments, which makes it highly flexible.
#' # If you don't like the default style of coefplot. No worries,
#' # you can set *your* default by using the function
#' # setFixest_coefplot()
#'
#' # Estimation
#' est = feols(Petal.Length ~ Petal.Width + Sepal.Length +
#' Sepal.Width | Species, iris)
#'
#' # Plot with default style
#' coefplot(est)
#'
#' # Now we permanently change some arguments
#' dict = c("Petal.Length"="Length (Petal)", "Petal.Width"="Width (Petal)",
#' "Sepal.Length"="Length (Sepal)", "Sepal.Width"="Width (Sepal)")
#'
#' setFixest_coefplot(ci.col = 2, pt.col = "darkblue", ci.lwd = 3,
#' pt.cex = 2, pt.pch = 15, ci.width = 0, dict = dict)
#'
#' # Tadaaa!
#' coefplot(est)
#'
#' # To reset to the default settings:
#' setFixest_coefplot("all", reset = TRUE)
#' coefplot(est)
#'
setFixest_coefplot = function(style, horiz = FALSE, dict = getFixest_dict(), keep,
ci.width = "1%", ci_level = 0.95, pt.pch = 20, pt.bg = NULL,
cex = 1, pt.cex = cex, col = 1:8, pt.col = col, ci.col = col,
lwd = 1, pt.lwd = lwd, ci.lwd = lwd, ci.lty = 1, grid = TRUE,
grid.par = list(lty = 3, col = "gray"), zero = TRUE,
zero.par = list(col = "black", lwd = 1), pt.join = FALSE,
pt.join.par = list(col = pt.col, lwd = lwd), ci.join = FALSE,
ci.join.par = list(lwd = lwd, col = col, lty = 2), ci.fill = FALSE,
ci.fill.par = list(col = "lightgray", alpha = 0.5),
ref.line = "auto",
ref.line.par = list(col = "black", lty = 2), lab.cex,
lab.min.cex = 0.85,
lab.max.mar = 0.25, lab.fit = "auto", xlim.add, ylim.add, sep, bg,
group = "auto", group.par = list(lwd = 2, line = 3, tcl = 0.75),
main = "Effect on __depvar__",
value.lab = "Estimate and __ci__ Conf. Int.",
ylab = NULL, xlab = NULL, sub = NULL, reset = FALSE){
fm_cp = formals(coefplot)
arg_list = names(fm_cp)
# arg_no_default = c("object", "se", "ci_low", "ci_high", "drop", "order", "ref", "add", "only.params", "as.multiple", "...", "x", "x.shift")
# m = fm_cp[!names(fm_cp) %in% arg_no_default]
# cat(gsub(" = ,", ",", paste0(names(m), " = ", sapply(m, deparse), collapse = ", ")))
iplot_default = list()
#
# Controls
#
check_arg(style, "character scalar")
check_arg(ci.width, "scalar(numeric, character) GE{0}")
check_arg(ci_level, "numeric scalar GT{0} LT{1}")
check_arg(lwd, ci.lwd, "numeric scalar GE{0}")
check_arg(grid, zero, "logical scalar")
check_set_arg("L0 list NULL{list()}", grid.par, zero.par, pt.join.par, ref.line.par)
check_arg(reset, "logical scalar")
#
# Code
#
if(missing(style)){
style = if(reset) "all" else "default"
}
if(style == "all" && reset){
opts = list(default = list(), iplot = iplot_default)
options("fixest_coefplot" = opts)
return(invisible(NULL))
}
opts = getOption("fixest_coefplot")
if(!is.list(opts)){
warning("Wrong format of getOption('fixest_coefplot'), the options of coefplot are reset.")
opts = list(default = list(), iplot = iplot_default)
}
if(reset){
if(style == "iplot"){
my_opt = iplot_default
} else {
my_opt = list()
}
} else {
my_opt = opts[[style]]
if(is.null(my_opt)){
my_opt = list()
} else if(!is.list(opts)){
warning("Wrong format of getOption('fixest_coefplot') for style '", style, "', the options of coefplot for this style are reset.")
my_opt = list()
}
}
mc = match.call()
all_args = setdiff(names(mc), c("", "reset"))
mc = mc[all_args]
# now we find which args for which we delay evaluation
arg_var_default = lapply(fm_cp[!names(fm_cp) %in% all_args], all.vars)
arg_var_default = arg_var_default[lengths(arg_var_default) > 0]
default.eval = sapply(arg_var_default, function(x) any(x %in% all_args))
if(any(default.eval)){
default.arg = names(default.eval)[default.eval]
mc[default.arg] = fm_cp[default.arg]
# we re order
my_fact = factor(names(mc), levels = arg_list)
mc = mc[order(my_fact)]
all_args = names(mc)
}
extra_values = my_opt$extra_values
if(!is.list(extra_values)) extra_values = list()
for(arg in all_args){
my_arg = mc[[arg]]
my_arg_vars = all.vars(my_arg)
if(length(my_arg_vars) == 0 || !(any(my_arg_vars %in% setdiff(arg_list, "dict")))){
if("par" %in% all.names(my_arg)){
my_opt[[arg]] = my_arg
} else {
my_opt[[arg]] = eval(my_arg)
}
} else {
my_extra_args = setdiff(my_arg_vars, arg_list)
if(length(my_extra_args) > 0){
for(v in my_extra_args) extra_values[[v]] = eval(str2lang(v), parent.frame())
}
# control that the order is proper
arg2eval = intersect(my_arg_vars, arg_list)
order_arg = which(arg_list == arg)
order_arg2eval = sapply(arg2eval, function(x) which(arg_list == x))
if(any(order_arg2eval > order_arg)){
qui = which(order_arg2eval > order_arg)[1]
stop("In argument ", arg, ": its value (", deparse(my_arg), ") is set with the argument ", arg2eval[qui], " which occurs after. If you set an argument with the value of another argument: you must use arguments appearing before.")
}
my_opt[[arg]] = my_arg
}
}
if(length(extra_values) > 0) my_opt$extra_values = extra_values
opts[[style]] = my_opt
options("fixest_coefplot" = opts)
}
#' @rdname setFixest_coefplot
getFixest_coefplot = function(){
getOption("fixest_coefplot")
}
####
#### plot method ####
####
#' plot methods for `fixest` and `fixest_multi` objects
#'
#' Plot method reporting the coefficient estimates and their confidence intervals.
#' This is a wrapper to the more complete functions [`coefplot`] and [`iplot`].
#'
#'
#' @inheritParams coefplot
#' @inheritSection etable Arguments keep, drop and order
#'
#' @param x A `fixest` estimation, for example from [`feols`].
#' @param ... Other arguments to be passed to [`coefplot`].
#'
#' @details
#' By default `plot.fixest` runs [`coefplot`] unless the estimation includes
#' the function [`sunab`], in which case it uses [`iplot`].
#'
#' The switch to `iplot` can be made with the argument `do_iplot = TRUE`.
#'
#' @return
#' It returns invisibly the data used to create the graph.
#'
#' @examples
#'
#' #
#' # Single estimation
#' #
#'
#' est = feols(Ozone ~ Temp + Solar.R, airquality)
#' plot(est)
#'
#' # focus only on the variables
#' plot(est, drop = "Cons")
#'
#' #
#' # Multiple estimations
#' #
#'
#' est_mult = feols(Ozone ~ csw(Temp, Solar.R, Wind), airquality)
#' plot(est_mult, drop = "Const")
#'
#' #
#' # DiD estimation: Sun & Abraham
#' #
#'
#' data(base_stagg)
#' # The DiD estimation
#' res_sunab = feols(y ~ x1 + sunab(year_treated, year) | id + year, base_stagg)
#' plot(res_sunab)
#'
#'
#'
#'
plot.fixest = function(x, vcov = NULL, add = FALSE, horiz = FALSE, do_iplot = NULL,
zero = TRUE, zero.par = TRUE, dict = NULL,
keep = NULL, drop = NULL, order = NULL,
ci.width = "1%", ci_level = 0.95, plot_prms = list(),
ylim = NULL, xlim = NULL,
pch = c(20, 17, 15, 21, 24, 22), col = 1:8, cex = 1, lty = 1, lwd = 1,
pt.pch = pch, pt.bg = NULL, pt.cex = cex, pt.col = col,
ci.col = col, pt.lwd = lwd, ci.lwd = lwd, ci.lty = lty,
main = "Effect on __depvar__",
value.lab = "Estimate and __ci__ Conf. Int.",
ylab = NULL, xlab = NULL, sub = NULL, ...){
check_arg(do_iplot, "NULL logical scalar")
if(is.null(do_iplot)){
do_iplot = FALSE
if(inherits(x, "fixest") && isTRUE(x$is_sunab)){
do_iplot = TRUE
} else if(inherits(x, "fixest_multi") && isTRUE(x[[1]]$is_sunab)){
do_iplot = TRUE
}
}
coefplot(objects = x, vcov = vcov, add = add, horiz = horiz, do_iplot = do_iplot,
zero = zero, zero.par = zero.par,
dict = dict, keep = keep, drop = drop, order = order, ci.width = ci.width,
ci_level = ci_level, plot_prms = plot_prms, ylim = ylim,
xlim = xlim, pch = pch, col = col, cex = cex, lty = lty, lwd = lwd,
pt.pch = pt.pch, pt.bg = pt.bg, pt.cex = pt.cex, pt.col = pt.col,
ci.col = ci.col, pt.lwd = pt.lwd, ci.lwd = ci.lwd, ci.lty = ci.lty,
main = main, value.lab = value.lab, ylab = ylab, xlab = xlab, sub = sub,
...)
}
#' @describeIn plot.fixest
plot.fixest_multi = plot.fixest
Try the fixest package in your browser
Any scripts or data that you put into this service are public.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Embedding an R snippet on your website
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.