R resources for GMM estimation and testing of multifactor asset pricing models

Question

Has anyone seen R script for GMM estimation and testing of asset pricing models such as Fama-French 3-factor or similar? Ideally, I would like to have R scripts corresponding to Cochrane "Asset Pricing" (2005) chapter 12 and related chapters.

Now, I have looked around and found some bits for the CAPM (including in the GMM package vignette) and some roughly equivalent bits elsewhere but not much else, and nothing for multifactor models. Perhaps R is not popular in the field of asset pricing, or perhaps I looked in the wrong places, I do not know. Any tips will be appreciated.

Answer 1

In this answer, I provide Matlab code for implementing a two-step GMM test of a multi-factor asset pricing model. I closely follow Cochrane (2005) which is an amazig book. See also this short video. The code begins from uploading the data and goes on to perform the GMM, calculate $t$-statistics and conduct the $J$-test. I hope it's easy enough to translate the code below to most other languages.

Upload the data

This first step is simple and just about uploading data from Ken French's website. For this example, I use the Fama-French five-factor model and the 25 portfolios sorted on size and book/market as test assets. My sample is July 1963 - June 2023.

clearvars                   % clear workspace
close all                   % close all windows
clc                         % clear command window

%%
% Upload raw data from Ken's website
TestAssets = readmatrix('25_Portfolios_5x5');                % 721 x 26 array
Factors = readmatrix('F-F_Research_Data_5_Factors_2x3');     % 721 x 7 array

% flip dimensions
TestAssets = TestAssets';
Factors = Factors';

% drop first row which would contain year/month
Factors(1,:) = [];          
TestAssets(1,:) = [];

% drop first column which would contain portfolio names
Factors(:,1) = [];          
TestAssets(:,1) = [];

% extract risk-free rate from factor file
Rf = Factors(end,:);
Factors(end,:) = [];

%Number of months, test assets and factors:
[N,T] = size(TestAssets);
[K,~] = size(Factors);

%Calculate excess returns
Rf = Rf/100;
Factors = Factors/100;
TestAssets = TestAssets/100 - Rf;

Note that $T=720$, $K=5$, and $N=25$. To check that everything looks alright, let's plot the average returns of our test assets

AvgRet = mean(TestAssets,2)*100;  % time series average 
AvgRet = reshape(AvgRet,[5,5]);
bar3(AvgRet)
xlabel('size')
ylabel('value')
zlabel('mean mo return %')

The result looks like this and indicates a positive value premium (amongst small stocks).

Do the GMM test

Let's turn to the fun bit. I implement the closed-form solution from Cochrane (2005) and follow his notation. In particular, see sections 10.1 and 13.2.

d = 1/T*TestAssets*Factors';
ERe = mean(TestAssets,2);
I = eye(N);  % identity matrix
% first step
    b_1 = inv(d'*d)*d'*ERe;
    SDF_1 = 1-b_1'*Factors;
    Errors_1 = SDF_1.*TestAssets;  % pricing errors of our model
    S_1 = cov(Errors_1'); % you can use fancier methods here (see section 11.7)
    W_1 = inv(S_1);   % weighting matrix for next step
    g_1 = mean(Errors_1,2); % average pricing errors (our moment conditions)
    cov_b_1 = 1/T*inv(d'*d)*d'*S_1*d*inv(d'*d);
    cov_g_1 = 1/T * (I-d*inv(d'*d)*d')*S_1*(I-d*inv(d'*d)*d');
    tstat_b_1 = b_1./diag(sqrt(cov_b_1));
% second step    
    b_2 = inv(d'*W_1*d)*d'*W_1*ERe;
    SDF_2 = 1-b_2'*Factors;
    Errors_2 = SDF_2.*TestAssets;
    S_2 = cov(Errors_2'); % you can use fancier methods here (see section 11.7)
    g_2 = mean(Errors_2,2);
    cov_b_2 = 1/T*inv(d'*inv(S_2)*d);
    cov_g_2 = 1/T * (S_2 - d*inv(d'*inv(S_2)*d)*d');
    tstat_b_2 = b_2./diag(sqrt(cov_b_2));
    J_2 = T*g_2'*inv(S_2)*g_2; %J-test statistic for overidentifying restrictions 
    J_crit_value = chi2inv(0.95,N-K); %H0: model is correct (all moments are zero)
    % reject H0 if J_2 > J_crit_value (ie the model doesn't price the test assets)

For example, the point estimate for SMB is 5.40 ($t$-statistic: 3.56). The $J$ test is rejected, however, suggesting that the model doesn't price the 25 portfolios (btw: the test is a really high hurdle and rarely passed).

Notes

$d$ is simply a second-moment matrix ($d = \mathbb{E}[R^e f] = \mathbb{C}\text{ov}(R^e,f)+\mathbb{E}[R^e]\mathbb{E}[f]$)

MyCov = NaN(N,K);
for i=1:N
    for j=1:K
        temp =  cov(TestAssets(i,:),Factors(j,:),1);
        MyCov(i,j) = temp(1,2)+mean(TestAssets(i,:))*mean(Factors(j,:));
    end
end
%compare d with MyCov!

We can implement a Do-it-Yourself GMM. Suppose we're not as clever as Cochrane (or Hansen) and do not know how to find the above closed-form solution. We can still use a naive brute-force approach to perform GMM. The aim is simply to minimise the pricing errors from our model. The objective function is $g(b)'Wg(b)$. The matlab code is

    W = eye(N);
    SDF = @(b) 1-b'*Factors;
    Errors = @(b) SDF(b).*TestAssets;
    g = @(b) mean(Errors(b),2);
    ObjFun = @(b) g(b)'*W* g(b);
    my_b_1 = fminsearch(ObjFun,ones(K,1));
    % compare my_b_1 with b_1 from the closed-form solution.

More notes. It's easy to compare different factor models using GMM (see Section 13.3 in Cochrane's textbook). Note that the standard errors from the above code are a tiny bit off because I do not correct for the fact that the factor mean is estimated from the same sample. Cochrane writes himself: In my experience so far with this method, the correction for the fact that $Ef$ is estimated is very small in practice, so that little damage is done in ignoring it (as is the case with the Shanken correction). So I kept it simple. He does explain how to fix the standard errors though (pages 257/258).