Using the ExPanDaR package for panel data exploration

Data Preparation

ExPanDaR expects the data to be organized in long format. This implies that each observation (a row) is identified by cross-sectional and time series identifiers and that variables are organized by columns. While you can have a vector of variables jointly determining the cross-section, the time-series needs to be identified by a unique sortable variable. The ExPanDaR functions treat cross-sectional identifiers as factors and expect the time-series identifier to be coercible into an ordered factor.

For this walk-through I will use the data set russell_3000, which comes with the package. It contains some financial reporting and stock return data of Russell 3000 firms and has been collected from Google Finance and Yahoo Finance using the tidyquant package in the summer of 2017. A word of caution: While the data appears to be relatively decent quality I would advise against using this data for scientific work without verifying its integrity first. These are the variables included in the data.

Let’s see whether coid and period identify a panel observation.

any(duplicated(russell_3000[,c("coid", "period")]))
## [1] FALSE

This seems to be the case.

Now you can start exploring the data interactively. The code

ExPanD(df = russell_3000, df_def = russell_3000_data_def,
       config_list = ExPanD_config_russell_3000,
       abstract = "Data is from Google Finance and Yahoo Finance.")

will start the ExPanD app (online instance can be found here).

But in this blog post I want to showcase the functions of the ExPanDaR package. So, as a first step, let’s use ExPanDaR’s function prepare_missing_values_graph() to eyeball how frequently observations are missing in the data set.

prepare_missing_values_graph(russell_3000, ts_id = "period")

OK. This does not look too bad. Only FY2013 seems odd, as some variables are completely missing. Guess why? They are calculated using lagged values of total assets. So, in the following, let’s focus on some profitability measures and on the fiscal years 2014 to 2016 (a short panel, I know). Time to check the descriptive statistics using the prepare_descriptive_table() function.

r3 <- droplevels(russell_3000[russell_3000$period > "FY2013",
                              c("coid", "coname", "period", "sector", "toas",
                                "nioa", "cfoa", "accoa", "return")])
t <- prepare_descriptive_table(r3)
t$kable_ret  %>%
  kable_styling("condensed", full_width = F, position = "center")

Take a look at the minima and the maxima of some of the variables (e.g., net income over assets (nioa)). Normally, it should be around -50 % to + 50%. Our measure has a minimum way below -50 %. One thing that comes very handy when dealing with outliers is a quick way to observe extreme values. prepare_ext_obs_table() might be helpful here.

t <- prepare_ext_obs_table(r3, cs_id = "coname", ts_id = "period", var = "nioa")
t$kable_ret %>%
  kable_styling("condensed", full_width = F, position = "center")

In a real life research situation, you might want to take a break and check your data as well as the actual financial statements to see what is going on. In most cases, you will see that the outliers are caused by very small denominators (average total assets in this case). To reduce the effect of these outliers on your analysis, you can winsorize (or truncate) them by using the treat_outliers() function.

r3win <- treat_outliers(r3, percentile = 0.01)
t <- prepare_ext_obs_table(r3win, cs_id = "coname", ts_id = "period", var = "nioa")
t$kable_ret %>%
  kable_styling("condensed", full_width = F, position = "center")

Descriptive Statistics

This looks better. Let’s look at the winsorized descriptive statistics.

t <- prepare_descriptive_table(r3win)
t$kable_ret  %>%
  kable_styling("condensed", full_width = F, position = "center")

Ok. Most firms have positive cash flows from operations and positive net income. Accruals (the difference between net income and cash flow from operations) are mostly negative. How does firms’ net income vary across sectors?

ret <- prepare_by_group_bar_graph(r3win, by_var = "sector", var = "nioa",
                                  stat_fun = mean, order_by_stat = TRUE)
## Warning: package 'bindrcpp' was built under R version 3.4.4
ret$plot

Whoa. This is not looking too good for the health sector. Let’s take a quick look at the distribution.

ret <- prepare_by_group_violin_graph(r3win, by_var = "sector", var = "nioa",
                                     order_by_mean = TRUE)
ret

You can see that there is a relatively large mass of loss making firms in the health sector. This is not too surprising given that many publicly-listed health firms are relatively small growth firms. Also, you can see the effect of regulated profits in the public utilities sector.

Time Trends

Additional visuals are available for exploring time trends. prepare_trend_graph() can be used for comparing variables…

graph <- prepare_trend_graph(r3win[c("period", "nioa", "cfoa", "accoa")], "period")
graph$plot

… and for eyeballing the distributional properties of a single variable over time you have prepare_quantile_trend_graph().

graph <- prepare_quantile_trend_graph(r3win[c("period", "return")], "period", c(0.05, 0.25, 0.5, 0.75, 0.95))
graph$plot

Nothing special going on here (not really surprising, given the short time span that the sample covers).

Correlation Analysis

I am sure that you won’t care but I am a big fan of correlation tables. prepare_correlation_table() prepares a table reporting Pearson correlations above and Spearman correlations below the diagonal.

t<- prepare_correlation_table(r3win, bold = 0.01, format="html")
t$kable_ret %>%
  kable_styling("condensed", full_width = F, position = "center")

In fact, I like correlations so much that especially for samples containing many variables I use prepare_correlation_graph() to display a graphic variant based on the corrplot package. See for yourself.

ret <- prepare_correlation_graph(r3win)

Scatter Plots

prepare_scatter_plot() produces the mother of all plots, the scatter plot. Here we take a look at the association of net income and stock returns.

prepare_scatter_plot(r3win, x="nioa", y="return", color="sector", size="toas", loess = 1)

They are mildly positively associated indicating that more profitable firms also tend to have higher stock returns. Do you see the structural break around nioa == 0? This indicates that firms and/or the stock market treats losses differently than gains. Researchers in the area of accounting tend to find that interesting.

Regression Tables

Finally, if you happen to be a fan of starred numbers, you can also quickly produce regression tables by using the function prepare_regression_table() that calls lfe::felm() for OLS and glm() for binary logit models. The tables are then constructed by calling stargazer::stargazer(), allowing for plain text, html and latex output.

You can construct tables by mixing different models…

dvs <- c("return", "return", "return", "return", "return", "return")
idvs <- list(c("nioa"), 
             c("cfoa"), 
             c("accoa"), 
             c("cfoa", "accoa"), 
             c("nioa", "accoa"), 
             c("nioa", "accoa")) 
feffects <- list("period", "period", "period", 
                 c("period", "coid"), c("period", "coid"), c("period", "coid"))
clusters <- list("", "", "", c("coid"), c("coid"), c("period", "coid"))
t <- prepare_regression_table(r3win, dvs, idvs, feffects, clusters)
htmltools::HTML(t$table)

… or by applying one model on different sub-samples.

t <- prepare_regression_table(r3win, "return", c("nioa", "accoa"), byvar="period")
htmltools::HTML(t$table)

Conclusion

This is all there is (currently). All these functions are rather simple wrappers around established R functions. They can easily be modified to fit your needs and taste. Take look at the github repository of the ExPanDaR package for the code and do not hesitate to touch base for feedback. Enjoy!