Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

implement TIS #22

Merged
merged 11 commits into from
Mar 17, 2024
Prev Previous commit
Next Next commit
changing effect default (with warning) and adding TIS tests (alongsid…
…e some testing changes)
  • Loading branch information
moritzpschwarz committed Mar 4, 2024
commit 62e61071b80fd2d0e87cb2a71024ddf24695346f
2 changes: 1 addition & 1 deletion R/identify_indicator_timings.R
Original file line number Diff line number Diff line change
Expand Up @@ -137,7 +137,7 @@ identify_indicator_timings <- function(object, uis_breaks = NULL){
output$csis <- csis
output$fesis <- fesis
output$cfesis <- cfesis
output$trends <- trends
output$tis <- trends
output$uis_breaks <- if(nrow(uis_indicators)>0) {uis_indicators} else{NULL}

return(output)
Expand Down
12 changes: 8 additions & 4 deletions R/isatpanel.R
Original file line number Diff line number Diff line change
Expand Up @@ -13,8 +13,8 @@
#' @param y Deprecated. The dependent variable. Can be used when data, index, and formula are not specified.
#' @param id Deprecated. Can be used when data, index, and formula are not specified. Must be a vector of the grouping variable as a character or factor
#' @param time Deprecated. Can be used when data, index, and formula are not specified. Must be a vector of the time variable as an integer or numeric.
#' @param mxreg The co-variates matrix
#' @param effect Fixed Effect specification. Possible arguments: "twoways", "individual", "time", or "none".
#' @param mxreg Deprecated.The co-variates matrix. Superseded by the formula argument.
#' @param effect Fixed Effect specification. Possible arguments: "twoways" (Default), "individual", "time", or "none".
#' @param na.remove remove NAs
#' @param engine Estimation function to use. Default is NULL, which uses the default estimation procedure of the gets package. Alternatives are "fixest", "plm", or "felm".
#' @param user.estimator Use a user.estimator
Expand All @@ -40,6 +40,7 @@
#' @param cfesis_var The cfesis method can be conducted for all variables (default) or just a subset of them. If you want to use a subset, please specify the column names of the variable in a character vector.
#' @param cfesis_id The cfesis method can be conducted for all individuals/units (default) or just a subset of them. If you want to use a subset, please specify the individuals/units to be tested in a character vector.
#' @param plot Logical. Should the final object be plotted? Default is TRUE. The output is a combination of \code{plot} and [plot_grid()] using the \code{cowplot} package.
#' @param print.searchinfo logical. If \code{TRUE} (default), then detailed information is printed.
#'
#' @return A list with class 'isatpanel'.
#' @export
Expand Down Expand Up @@ -87,7 +88,7 @@
data=NULL,
formula=NULL,
index=NULL,
effect = c("individual"),
effect = c("twoways"),

na.remove = TRUE,
engine = NULL,
Expand All @@ -113,6 +114,7 @@
t.pval = 0.001,

plot = TRUE,
print.searchinfo = TRUE,
plm_model = "within",

y=NULL,
Expand All @@ -125,6 +127,7 @@

# Error checks
if (!effect %in% c("twoways", "individual", "time","none")) {stop("Error in Fixed Effect Specification (effect). Possible values for effect are: 'twoways', 'individual', 'time', or 'none'.")}
if (missing(effect) & print.searchinfo){warning("New default for effect in 'isatpanel': Used to be 'individual', now 'twoways'. To quiet this message provide the argument 'effect' or select 'print.searchinfo = FALSE'.")}

if ((effect == "both" | effect == "time") & jiis == TRUE) {stop("You cannot use time fixed effects and jiis = TRUE. These would be perfectly collinear. Either set jiis = FALSE or use effect = 'individual'.")}

Expand Down Expand Up @@ -473,7 +476,7 @@
# remove any duplicate columns
drop <- union(which(duplicated(as.list(current),fromLast = TRUE)),which(duplicated(as.list(current),fromLast = FALSE)))
if (!identical(drop,integer(0))) {
current <- current[,-drop]

Check warning on line 479 in R/isatpanel.R

View check run for this annotation

Codecov / codecov/patch

R/isatpanel.R#L479

Added line #L479 was not covered by tests
}

# Add to Breaklist (uis list)
Expand Down Expand Up @@ -598,7 +601,8 @@

options(mc.warning = FALSE)
#sis <- FALSE # don't allow sis argument - does not make sense in a panel context, only JSIS makes sense
ispan <- gets::isat(y, mxreg = mx, iis = iis, sis = FALSE, uis = sispanx, user.estimator = user.estimator, mc = FALSE, t.pval = t.pval, ...)
ispan <- gets::isat(y, mxreg = mx, iis = iis, sis = FALSE, uis = sispanx, user.estimator = user.estimator, mc = FALSE, t.pval = t.pval,
print.searchinfo = print.searchinfo, ...)
#ispan <- isat.short(y, mxreg = mx, iis=iis, sis=FALSE, uis=sispanx, user.estimator = user.estimator, mc=FALSE, ...)

###############################
Expand Down
4 changes: 2 additions & 2 deletions R/plot.isatpanel.R
Original file line number Diff line number Diff line change
Expand Up @@ -58,8 +58,8 @@ plot.isatpanel <- function(x, max.id.facet = 16, facet.scales = "free", title =
}

# tis
if(!is.null(df_identified$trends)){
g = g + geom_vline(data = df_identified$trends, aes(xintercept = .data$time,color="lightblue"))
if(!is.null(df_identified$tis)){
g = g + geom_vline(data = df_identified$tis, aes(xintercept = .data$time,color="lightblue"))
}

# cfesis
Expand Down
9 changes: 6 additions & 3 deletions man/isatpanel.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

21 changes: 16 additions & 5 deletions tests/testthat/test-1-isatpanel-short.R
Original file line number Diff line number Diff line change
Expand Up @@ -23,12 +23,19 @@ pandata_simulated <- pandata_simulated[pandata_simulated$year>1979,]
# Unit testing
test_that("Initial Tests Isatpanel on simulated data",{

expect_message(isatpanel(data = pandata_simulated,formula = gdp~temp, index = c("country","year"),fesis = TRUE))
expect_message(isatpanel(data = pandata_simulated,formula = gdp~temp, index = c("country","year"), fesis = TRUE, effect = "twoways",
print.searchinfo = TRUE))

#newmethod <- isatpanel(data = pandata_simulated,formula = gdp~temp + I(temp^2), index = c("country","year"),fesis=TRUE)

expect_message(isatpanel(data = pandata_simulated,formula = gdp~temp + I(temp^2),
index = c("country","year"),fesis=TRUE, ar = 1))
index = c("country","year"),fesis=TRUE, ar = 1, effect = "twoways"))


expect_warning(isatpanel(data = pandata_simulated,formula = gdp~temp + I(temp^2),
index = c("country","year"),fesis=TRUE, ar = 1),
regexp = "New default for effect in 'isatpanel': Used to be 'individual', now 'twoways'.")


})

Expand All @@ -50,7 +57,9 @@ test_that("Isatpanel Test that missing values are removed",{
test_that("Test the cfesis and csis arguments",{

expect_silent(isatpanel(data = pandata_simulated,formula = gdp~temp + I(temp^2), index = c("country","year"),fesis=TRUE, cfesis = TRUE, ar = 1, print.searchinfo = FALSE))
expect_silent(newmethod_cfesis_sub <- isatpanel(data = pandata_simulated,formula = gdp~temp + I(temp^2), index = c("country","year"),fesis=TRUE, cfesis = TRUE,cfesis_id = c("2","3"), ar = 1, print.searchinfo = FALSE))
expect_silent(newmethod_cfesis_sub <- isatpanel(data = pandata_simulated,formula = gdp~temp + I(temp^2),
index = c("country","year"),fesis=TRUE, cfesis = TRUE,cfesis_id = c("2","3"), ar = 1,
print.searchinfo = FALSE))

})

Expand All @@ -62,8 +71,10 @@ test_that("Test the cfesis and csis arguments",{

test_that("Standard Error Options using fixest",{
skip_on_cran()
expect_silent(result <- isatpanel(data = pandata_simulated,formula = gdp~temp + I(temp^2), index = c("country","year"),fesis=TRUE, ar = 1, print.searchinfo=FALSE,engine = "fixest"))
expect_silent(result <- isatpanel(data = pandata_simulated,formula = gdp~temp + I(temp^2), index = c("country","year"),fesis=TRUE, ar = 1, print.searchinfo=FALSE,engine = "fixest", cluster = "individual"))
expect_silent(result <- isatpanel(data = pandata_simulated,formula = gdp~temp + I(temp^2), index = c("country","year"),fesis=TRUE, ar = 1,
print.searchinfo=FALSE,engine = "fixest"))
expect_silent(result <- isatpanel(data = pandata_simulated,formula = gdp~temp + I(temp^2), index = c("country","year"),fesis=TRUE, ar = 1,
print.searchinfo=FALSE,engine = "fixest", cluster = "individual"))
})


Expand Down
File renamed without changes.
122 changes: 122 additions & 0 deletions tests/testthat/test-5-tis.R
Original file line number Diff line number Diff line change
@@ -0,0 +1,122 @@

# # Single TS Example
# # Generate some random data
# set.seed(123)
# x <- rnorm(50, mean = 100)
# ep <- rnorm(50, sd = 0.2)
# trend <- 1951:2000
# trendbreak <- c(rep(0,24),1:26) # impose a trendbreak from 1975
# y <- 10 + 0.5 * x + 0.1 * trend - 0.2 * trendbreak + ep
# df <- data.frame(id = "A",
# year = 1951:2000,
# y = y,
# x = x,
# trend = trend,
# trendbreak = trendbreak,
# ep = ep)
#
# # Show a model without considering the trendbreak
# #gets::arx(y = y, mc = TRUE, mxreg = df[,c("x","year")], plot = TRUE)
#
# # Running TIS
# #gets::isat(y = y, mc = TRUE, mxreg = df[,c("x","year")], sis = FALSE, tis = TRUE, plot = TRUE)
#
#
# # show this result as well using two break indicators
# # this first one always starts 0,0,1,2...
# num_trend <- trend - 1950
# num_trend_break <- rep(0,50)
# num_trend_break[25:50] <- 1:26
#
# # this one works just like MIS or csis/cfesis in getspanel
# # the original trend would have been 1,2,3,4
# # the MIS result for this would be 0,0,3,4
# full_trend_break <- 1:50
# full_trend_break[1:24] <- 0
#
# df$full_trend_break <- full_trend_break
# df$num_trend_break <- num_trend_break
#
# # the MIS approach obviously does not produce the right result
# #gets::arx(y = y, mc = TRUE, mxreg = df[,c("x","year","full_trend_break")], plot = TRUE)
#
# # while the TIS approach works
# #gets::arx(y = y, mc = TRUE, mxreg = df[,c("x","year","num_trend_break")], plot = TRUE)

###
# Generate Three Unit Example
###
set.seed(123)
# Generate some random data for the two control countries
xA <- rnorm(50, mean = 100)
xB <- rnorm(50, mean = 30)
xC <- rnorm(50, mean = 70)

epA <- rnorm(50, sd = 0.2)
epB <- rnorm(50, sd = 0.2)
epC <- rnorm(50, sd = 0.2)

trend <- 1951:2000
trendbreak <- c(rep(0,19),1:31) # impose a trendbreak from 1975

yA <- 10 + 0.5 * xA + 0.2 * trend - 0.3 * trendbreak + epA
yB <- 0.5 * xB + 0.1 * trend + epB
yC <- 0.5 * xC + 0.1 + epC

trial_df <- data.frame(year = rep(1951:2000,3),
id = c(rep("A",50),rep("B",50),rep("C",50)),
x = c(xA, xB, xC),
y = c(yA,yB,yC))


# Introduce a step shift in A from 40
trial_df_step <- trial_df
trial_df_step$y[40:50] <- trial_df_step$y[40:50]*1.02


test_that("TIS works", {
expect_silent(m1 <- isatpanel(trial_df, formula = y ~ x, index = c("id","year"), tis = TRUE, plot = FALSE, print.searchinfo = FALSE))
expect_type(get_indicators(m1), type = "list")
expect_identical(names(get_indicators(m1)), "tis")
expect_identical(get_indicators(m1)$tis$name, c("tisA.1952", "tisA.1970", "tisC.1956"))



expect_silent(m2 <- isatpanel(trial_df_step, formula = y ~ x, index = c("id","year"), effect = "twoways", tis = TRUE, print.searchinfo = FALSE)) # TIS approximates step shift
expect_silent(m3a <- isatpanel(trial_df_step, formula = y ~ x, index = c("id","year"), effect = "twoways", fesis = TRUE, print.searchinfo = FALSE)) # Time FE approximates trend (causes problems in C)
expect_silent(m3b <- isatpanel(trial_df_step, formula = y ~ x, index = c("id","year"), effect = "individual", fesis = TRUE, print.searchinfo = FALSE)) # Step Shift approximates trend (esp in B)
expect_silent(m4 <- isatpanel(trial_df_step, formula = y ~ x, index = c("id","year"), effect = "twoways", fesis = TRUE, tis = TRUE, print.searchinfo = FALSE)) # correct specification

expect_identical(get_indicators(m4)$fesis$name, "fesisA.1990")
expect_identical(get_indicators(m4)$tis$name, c("tisA.1955", "tisA.1970","tisB.1962"))



# this is a known issue that we are looking into
# here, the IIS are retained despite the fesis being active
# the result is valid, but not very efficient (reduces power)
# isatpanel(trial_df_step, formula = y ~ x, index = c("id","year"), fesis = TRUE, tis = TRUE, iis = TRUE)




# make sure tis_id works
expect_silent(m4a <- isatpanel(trial_df_step, formula = y ~ x, index = c("id","year"), effect = "individual",
fesis = TRUE, tis = TRUE, print.searchinfo = FALSE,
tis_id = c("A","C")))
expect_identical(get_indicators(m4a)$fesis$name, c("fesisA.1990","fesisB.1962","fesisB.1981"))
expect_identical(get_indicators(m4a)$tis$name, c("tisA.1952", "tisA.1975"))


expect_silent(m4b <- isatpanel(trial_df_step, formula = y ~ x, index = c("id","year"), effect = "individual",
fesis = TRUE, tis = TRUE, print.searchinfo = FALSE,
tis_id = c("A","C"),
fesis_id = c("A","C")))

expect_identical(get_indicators(m4b)$fesis$name, "fesisA.1990")
expect_identical(get_indicators(m4b)$tis$name, c("tisA.1952", "tisA.1975"))
})




Loading