Posted Question for R Users

I recently undertook a project where a colleague had about 12 .csv files that they wanted to merge. Each file had a common (key) variable 'Partner' (which is trading partner) with differing columns (variables) except for the common key variable. Actually, the first column (Partner) was the same in each data set, with each data set having 211 rows. The remaining columns were unique to each data set. I required a combined data set.

In SAS (PROC SQL) this would be a matter of a series of left joins. Working in the R environment, I executed what was equivalent to left joins via the 'merge' function. However this was tedious, only being able to join 2 tables at a time. In SAS I could use the merge function, which would allow me to merge all 12 tables in 1 data step.

Does anyone know of a better way to do this in R, as opposed to my merge statements? (see sample code below)

R Merge Statements:

# -------------------------------------------------------------------
#|  merge data sets with R merge function              
# ------------------------------------------------------------------
 
#  left join b_dat onto a_dat on variable Partner
 
 ab <- merge(a_dat,b_dat, by.a_dat = Partner, by.b_dat =Partner, all = FALSE, all.x = TRUE, all.y = FALSE)
 
# left join c_dat onto ab on variable Partner
 
 a_c <- merge(ab,c_dat, by.ab = Partner, by.c_dat =Partner, all = FALSE, all.x = TRUE, all.y = FALSE)
 
# I have about 10 more data sets to left join with a_c, is there a better way to join these 
# in R as opposed to pairwise merges like above?

Created by Pretty R at inside-R.org

I found some help on r-wiki, it merges all 3 data sets at once, but gives me extra redundant columns with .x and .y appended to their names. I'm not sure about these results.

my.list <- list(a_dat, b_dat, c_dat)
 
DF <- a_dat
for ( .df in my.list ) {
  DF <-merge(DF,.df,by.x="Partner", by.y="Partner", all = FALSE, all.x = TRUE, all.y = FALSE)
 }
 
dim(DF)

Created by Pretty R at inside-R.org

SAS - similar code for 3 data sets DAT_A, DAT_B, DAT_C

DATA ALL_COUNTRIES;
MERGE DAT_A (IN=A) DAT_B (IN=B)DAT_C (IN=C);
BY PARTNER;
IF B AND C; /* Actually its been so long since */
RUN; /* I have used a data step vs PROC SQL I'm not sure this*/
/* statement gives me a left join, but I'm not sure I can */
/* do any better than 2 at a time left joins in PROC SQL */
/* so that would not be any better than my R code above */

Data Science

 (Venn Diagram by Drew Conway)

Below is a great article on the rise of data science. Given my background in applied economics, and my current work experience in hacking and data mining, I find this article very interesting and accurate.

Article: 'The Rise of Data Science' .

"Data scientists don't have to be hard-core computer geeks, says Conway, but they do have to know their way around the IT landscape because that's where the data lives. Hacking abilities are important because data tends to reside in multiple locations, and in multiple systems. Finding and retrieving data sometimes requires the skills of a burglar -- even when the data is in the public domain, owned by your organization, or owned by another organization that has agreed to let you use it."

"The trend is toward a multi-disciplinary approach to extracting value from data. It's not just about math anymore. You also need technology skills, but what ultimately separates the analyst from the scientist is the dimension of artistic creativity. It's the soft skills that make the big difference."

The Venn diagram "is useful, but I don't think that one person can do all of this," says Elashoff. An interdisciplinary team, however, could possess the skills depicted in the diagram. "It's pretty hard to be an expert in all these areas."

Forecasting and Time Series Analysis

Time Series Analysis

Exponential Smoothing:

S = α Y _t +  (1-α) Y_t-1              
S =  portion of current average + portion of previous average = smooth value

Classical Equation:

Forecasted Value = T x S x C                        
FV =trend*seasonal index*cyclical index} 

Regression Equation:

T + S + C    
Y = trend regression with seasonal and cyclical dummies

Seasonal Indices: 

AVG_{season i}  /  AVG_total  = seasonal index for season i.

Econometric Forecasting Models

Distributed Lag Models

Y_t = a₀ + b₀ X _t +  λY_t-1  +  u_t

where λ is on the interval (0,1)  
 

e_t = r e_t-1  +  u_t    (serially correlated error terms)

Serially correlated error terms can be corrected using GLS, iterated maximum likelihood, instrumental variables.

Non-stationary: upward or downward trend in data, this can be corrected by take first differences:

           
         Y* = Y­_t – Y_t-1

Application: To predict  Y_t+1 use the current value of Y.  Predicting Y_t+2 can be predicted as a conditional forecast, given a forecast of Y_t+1 

ARIMA Models

Combines moving average (MA) and autoregressive (AR) models as well as differencing.

ARIMA(p,d,q)  where p = # of autoregressive terms
                                    d = # of 1st differences
                                    q = # moving average terms

Application:  ARIMA(2,1,1) :    Y*_t = b₀ + b₁Y*_t-1  + b₂Y*_t-2  + b₃ e_t-1  + e_t

2 autoregressive terms,  1 1^st difference, 1 moving average term

Seasonal effects can also be incorporated:

ARIMA(p,d,q)(P,D,Q)


where P = # of seasonally differenced autoregressive terms
          D = # of seasonal differences
          Q = # of seasonally differenced moving average terms


This model can also be represented as such:

[θ(B) θ(B^s) / Φ(B) Φ(B^s) Δ^dΔ^s  ] *   a_t 

where  

B = backshift operater such that BY = Yt-1

θ(B) = MA(q)

θ(B^s) = seasonal MA(Q)

Φ(B) = AR(p)

Φ(B^s) = seasonal AR(P)

Δ^d=  differencing

Δ^s = seasonal differencing

a_t = white noise term

References:

Using Econometrics: A Practical Guide (4th Edition)

A.H. Studenmund

Analysis of Variance

BIOMETRICS

COMPARING TWO MEANS: Xbar₁ – Xbar₂ / S_x1 – S_x2 ~ t 

can’t define numerator for >2 means

ANALYSIS OF VARIANCE: Comparing more than two means. When μ₁ not = μ₂, an estimate of variance using means will include a contribution attributable to the difference between population means. This ‘treatment’ effect will lead to a biased estimate of variance and large F.

F = σ ²+ trt / σ = S²_b  / S²_w ~ F(.)

If means are different, then  S²_b  > S²_w  , then F will be large. A large significant F value allows you to reject the null hypothesis that the means are equal. 

Linear Additive Model: Y_ij = U + T_i + E_ij

Or in matrix form: Y = Xb + e  where  Y = a matrix of individual observations. X = a ‘design matrix of 1’s and 0’s, & b =  a matrix of co-efficients of the population mean and treatment effects ‘t’.  or ( U + T_i) = Xb

SS Total = Y’Y – Y’CY’C(1/n)

SS Between/Model=b’X’Y - Y’CY’C(1/n)

SS Within/Error  =Y’Y’ – B’X’Y

Where ‘C’ is a matrix of 1’s

RCB:

                                    ( Blocks)

(treatments)               A            B            C             D             E

1            x             x            x            x            x

2            x            x            x            x            x

3            x            x            x            x            x

4            x            x            x            x            x

1.  If treatments or blocks differ, look at individual differences in means via LSD.
2.  Contrasts: Q = CiYi  ( orthogonal)   MS(Q) = Q² / r sum(Ci²), a subset of a linear function, a way of comparing different means.

LATIN SQUARE:

A	D	C	B
B	C	A	D
D	A	B	C
C	B	D	A

Each row is a comparable block with treatments A-D occurring only once in each block.

FACTORIAL EXPERIMENTS:

1. # treatments consists of several categories, levels
2. May be used to determine optimal levels
3. Interaction: difference in response levels-non additive,  measures homogeneity of  a response.

Polynomial Responses:  partitioning SS into linear and quadratic components- response surfaces

SPLIT PLOT DESIGNS

Whole plots divided into sub-plots where levels of factors are applied.

Block #1

A4B2 A4B1

A1B2 A1B1

A2B1 A2B2

A3B1 A3B2

A4B2 A4B1

Incomplete block re ‘A’ , complete block re ‘B’

SPLIT BLOCK EXPERIMENTS

Split plot in time vs. space ex:  time or cutting represents an abstract split plot ‘a’ cultivars, ‘b’ cuttings.

ANALYSIS OF COVARIANCE

1. Observed variation in Y is partly due to variation in X.
2. Use regression to eliminate effects that cannot be controlled by experimental design.

Ex: animals in a block with varying initial weights, use covariance to remove effect from experimental error.

Y_ij = μ + T_i + p_j + B(X_ij – Xbar) + e_ij

B = Cov(X,Y) / Var(X)

Carry out AOV on values adjusted for regression on an independent variable.

References:

Principles and Procedures of Statistics: A Biometrical Approach 3rd Edition. Robert George Douglas Steel, James Hiram Torrie, David A. Dickey

Mixed, Fixed, and Random Effects Models

The General Mixed Model:

Y  = Xβ + Z α+ e            (1)

β =  fixed effects co-efficient vector

X = fixed effects model matrix

α = random effects co-efficient matrix

Z =  random effects model matrix

The fixed part of the model is specified by Xβ and the random part by Zα+ e.

 The ‘random effects’ matrix (α) represents random effects that vary across individuals vs. the ‘fixed effects’ matrix (β) that represents effects that are the same across all individuals.

Panel Data: Cross sectional time series data, in most cases looking at hundreds or thousands of individuals (units) observed at several points across time, i.e. multiple observations per unit across time.

Heterogeneity Effect: Whether or not effects, or responses of individuals are the same across time, or if there are group differences.  If effects are not the same, and they are not accounted for, estimation errors result. Fixed and random effects models attempt to capture the heterogeneity effect. 

Given

y_it =b x_it + α_i + u_it                  (4)

 α_i is an unobserved individual effect

Fixed Effects Model (FE): α_i is correlated with x. In this model individual effects, or differences across individuals can be captured by shifts in the regression equation, or dummy variables.

As such we can estimate the fixed effects model as a Least Squares Dummy Variable model (LSDV):

y= Xb + dα + e            (5)

where d is a vector of dummy variables for each individual or unit effect.

Random Effects Model (RE): α_i is uncorrelated with x. Individual effects are randomly distributed across units.

y_it =b x_it + α_i+ u_it            (6)

Note that (4) and (6) are identical, the only difference is the assumptions about the individual effect α_i and how the models are estimated.

Note also that (6) differs from (1) in that we don’t have the matrix Z. The Random Effects model is in fact a special case of the general mixed model with a random intercept (α_i).

The random effects model is estimated using Generalized Least Squares (GLS) :

β_GLS = (X’Ω^-1X)^-1(X’Ω^-1Y) where Ω = I ⊗ Σ    (7)

Where Σ is the variance α_i+ u_it . If  Σ is unknown, it is estimated, producing a feasible generalized least squares estimate β_FGLS

Whether a FE or RE model should be used can be determined based on the Housman test (see Greene, for more details)

References:

Greene, Econometric Analysis, 5^th Edition

Instrumental Variables (IVs)

In a previous post I discussed how IVs and 2SLS can be used to correct for simultaneity bias in systems of simultaneous equations. 2SLS and IVs may also be useful in other contexts, such as the case of measurement error. Greene gives the following example:

Suppose you are estimating a consumption function with income as an explanatory variable:

C = f(I) + e

In general, it may be that income is not accurately reported. As a result, the estimation of C will suffer from measurement error. Greene, in the example given in his text, proposes that instead of using income as an explanatory variable, if we had data on the number of checks written per household, this variable would likely be correlated with income, but uncorrelated with measurement error. The number of checks written would be an instrumental variable.

IV Estimation via 2SLS:

In the context of the example above, let C be our dependent variable and x be income, then our regression equation becomes:

y = b x + e

An instrumental variable  ‘z’ is one that is correlated with ‘x’ but not ‘e’

E(y|z) = bE(x|z) + E(e|z), and by assumption E(e|z) =0

Stage 1: Regress x on z to get x*

Stage 2: Regress y on x*

This gives us b_IV = (z’x)^-1z’y and the result  y = b_IV x* +e

An important concept in IVs is the exclusion principle which states that the only way z impacts y is through z’s effect on x. In other words, the causal model is:

z→ x→ y

Related to this, I believe, is a statement you will find in 'Mostly Harmless Econometrics' by Angrist and Pischke (where s is the independent variable being instrumented):

"Intuitively, conditional on covariates, 2SLS retains only the variation in s that is generated by quasi-experimental variation- that is generated by the instrument z"

References:

Angrist and Pischke, Mostly Harmless Econometrics, 2009 

Greene, Econometric Analysis. 5^th Edition

See also this great blog post from Dr. Andrew Gelman with comments from Hal Varian: How to think about instrumental variables when you get confused

“Suppose z is your instrument, T is your treatment, and y is your outcome. So the causal model is z -> T -> y……. when I get stuck, I find it extremely helpful to go back and see what I've learned from separately thinking about the correlation of z with T, and the correlation of z with y. Since that's ultimately what instrumental variables analysis is doing.”

Hal Varian:

"You have to assume that the only way that z affects Y is through the treatment, T. So the IV model is
T = az + e
y = bT + d

It follows that
E(y|z) = b E(T|z) + E(d|z)

Now if we
1) assume E(d|z) = 0
2) verify that E(T|z) != 0
we can solve for b by division. Of course, assumption 1 is untestable.

An extreme case is a purely randomized experiment, where e=0 and z is a coin flip."

Simultaneous Equations, 2SLS, & IVs

Simultaneous Equations

(see also Instrumental Variables)

Simplifying the notation, let’s assume we have the following system of simultaneous/structural equations:

Y₁ = f(x₁,Y₂) +e₁  (1)

Y₂ = f(x₂,Y₁) +e₂  (2)

Both equations are necessary to characterize the relationships between the variables, they must be viewed as an entire system to grasp all of the feedback mechanisms involved. 

endogenous variables: those jointly determined by the system, in this case we have two endogenous variables, Y₁ and Y₂.

exogenous variables: those variables that are not jointly determined by the system, in this case x₁ and x₂. In othe words these variables are determined outside the system, or ‘exogenous’ to the system.

predetermined variables: exogenous variables and lagged endogenous variables (there are not lagged endogenous variables in this system, but an example would be Y_1,t-1.)

simultaneity bias: In the system represented above, we have the following consequence- a change in e₁ leads to a change in Y₁, but because Y₁ is a predictor of Y₂, this causes a change in Y₂, which feeds back into equation (1) causing a change in Y₁.  We end up with the following result:

Cov(Y₂,e₁) ≠  0

Recall the assumptions from classical regression require that the correlation between the error term and explanatory variables be zero. The condition of simultaneity bias above violates that assumption.  The estimated regression co-efficients become biased:

E(b) ≠  b

This can be corrected using 2-stage least squares (2SLS) and instrumental variables.

Two-Stage Least Squares

If we could find a variable ‘z’ highly correlated with Y₂ but uncorrelated with the error term e₂ we could run regressions on (1) and (2) and avoid the problem of simultaneity bias.  This variable ‘z’ is referred to as an ‘instrumental variable.’

Studenmund defines 2SLS as ‘a method of systematically creating instrumental variables to replace endogenous variables where they appear as explanatory variables in simultaneous equation systems.’

Equations (1) and (2) can be re-written in what’s referred to as reduced form:

Y₁ = f(x₁,x₂) +e₁  (3)

Y₂ = f(x₂,x₁) +e₂  (4)

Stage 1: Run OLS on (3) and (4) to get estimates  z₁ and z₂ for Y₁  and Y₂ respectively. 

z₁ = f(x₁,x₂) +e₁  (5)

z₂ = f(x₂,x₁) +e₂  (6)

Stage 2: Run OLS on (1) and (2) using instruments z₁ and z₂ to replace Y₁ and Y₂ where they appear as explanatory variables.

Y₁ = f(x₁,z₂) +e₁  (7)

Y₂ = f(x₂,z₁) +e₂  (8)

References:

Using Econometrics, Studenmund  (2001)