Processing math: 100%

In my earlier post from March 2018, I did introduce the rollapply function that executes a function on a rolling window basis. While this function is very useful, it needs a little modification for users to apply other general operations.

Originally, I faced this issue when I tried to compute the correlation matrix across different asset returns on a rolling window. For the demonstration, let’s consider the returns for all sector ETFs excluding real estate:

library(quantmod)
v <- c("XLE","XLU","XLK","XLB","XLP","XLY","XLI","XLV","XLF")
t1 <- "1990-01-01"
P.list <- lapply(v,function(x) get(getSymbols(x,from = t1)) )
P.list <- lapply(P.list,function(x) x[,6])
P <- Reduce(merge,P.list)
names(P) <- v
R <- na.omit(P/lag(P) - 1)

By default, rollapply executes the given function on each time series separately and returns a time series object. For instance,

tail(rollapply(R,25,mean))
                   XLE           XLU          XLK         XLB           XLP         XLY         XLI          XLV         XLF
2019-01-23 0.000405934 -1.762624e-03 0.0003577608 0.001246882 -0.0009785855 0.001524203 0.001472864 0.0005724254 0.002758011
2019-01-24 0.001347467 -2.660886e-04 0.0015414921 0.001719077 -0.0005482073 0.002707476 0.002332465 0.0010333171 0.003247331
2019-01-25 0.002800876 -5.444870e-04 0.0017911008 0.002256668 -0.0002366034 0.002836724 0.002641057 0.0014315231 0.003723034
2019-01-28 0.002870862 -6.661247e-04 0.0020484144 0.002740595  0.0003212731 0.003396979 0.003012958 0.0015001508 0.004101841
2019-01-29 0.004101322 -5.310597e-04 0.0023354609 0.003691053  0.0011276797 0.004034883 0.004298240 0.0022585259 0.004318732
2019-01-30 0.005052890 -5.192435e-05 0.0047939592 0.004476206  0.0018114857 0.005728323 0.005678859 0.0033033039 0.005350848

returns the 25 moving average for each one separately. On the other hand, if I try to compute the moving correlation, instead, I get the following

tail(rollapply(R,25,cor))
           XLE XLU XLK XLB XLP XLY XLI XLV XLF
2019-01-23   1   1   1   1   1   1   1   1   1
2019-01-24   1   1   1   1   1   1   1   1   1
2019-01-25   1   1   1   1   1   1   1   1   1
2019-01-28   1   1   1   1   1   1   1   1   1
2019-01-29   1   1   1   1   1   1   1   1   1
2019-01-30   1   1   1   1   1   1   1   1   1

which computes the correlation with the same ETF rather than other ETFs - as it treats each time series separately. As a remedy, one should add by.column = F argument to the rollapply function. In this case, the function returns a time series xts object, however, with 9×9=81 columns, where each column corresponds to the pairwise correlation between the 9 sector ETFs rather than a squared matrix.

COR <- rollapply(R,25,cor,by.column = F)
dim(COR)
[1] 5057   81
class(COR)
[1] "xts" "zoo"

What left to be done is to stack these vectors back into a correlation matrix, one for each time period. To do so, I will refer to the plyr package. The plyr package allows users to take an array (a), a data frame (d), or a list (l), execute a given function over the given object, and output the results in either format. For our case, I will input the time series COR object as an array and output it as a list, where each element in the list corresponds to the moving correlation matrix.

library(plyr)
COR.list <- alply(COR,1,function(x) matrix(x,nrow = ncol(R),  byrow = T ))

The second argument in the alply specifies the margin, where 1 indicates that the given function to be executed over the rows, while 2 states that it should be executed over the columns instead. The third argument, which takes a function, stacks each row of the COR object into a squared matrix. As a result, we have:

round(COR.list[[25]],2)
      [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9]
 [1,] 1.00 0.29 0.41 0.52 0.41 0.48 0.57 0.42 0.69
 [2,] 0.29 1.00 0.30 0.05 0.34 0.52 0.11 0.40 0.23
 [3,] 0.41 0.30 1.00 0.32 0.50 0.60 0.45 0.71 0.66
 [4,] 0.52 0.05 0.32 1.00 0.31 0.52 0.80 0.45 0.59
 [5,] 0.41 0.34 0.50 0.31 1.00 0.71 0.52 0.50 0.57
 [6,] 0.48 0.52 0.60 0.52 0.71 1.00 0.48 0.54 0.74
 [7,] 0.57 0.11 0.45 0.80 0.52 0.48 1.00 0.56 0.61
 [8,] 0.42 0.40 0.71 0.45 0.50 0.54 0.56 1.00 0.58
 [9,] 0.69 0.23 0.66 0.59 0.57 0.74 0.61 0.58 1.00

which is identical to correlation matrix computed over the first 25 days in the data

round(cor(R[1:25,]),2)
     XLE  XLU  XLK  XLB  XLP  XLY  XLI  XLV  XLF
XLE 1.00 0.29 0.41 0.52 0.41 0.48 0.57 0.42 0.69
XLU 0.29 1.00 0.30 0.05 0.34 0.52 0.11 0.40 0.23
XLK 0.41 0.30 1.00 0.32 0.50 0.60 0.45 0.71 0.66
XLB 0.52 0.05 0.32 1.00 0.31 0.52 0.80 0.45 0.59
XLP 0.41 0.34 0.50 0.31 1.00 0.71 0.52 0.50 0.57
XLY 0.48 0.52 0.60 0.52 0.71 1.00 0.48 0.54 0.74
XLI 0.57 0.11 0.45 0.80 0.52 0.48 1.00 0.56 0.61
XLV 0.42 0.40 0.71 0.45 0.50 0.54 0.56 1.00 0.58
XLF 0.69 0.23 0.66 0.59 0.57 0.74 0.61 0.58 1.00

Finally, one can either keep the rolling correlation matrix in a list or transform it back a time series using certain computation, e.g. construct portfolio weights and compute the out-of-sample return as a time series. As a finall demostration, I will show how one can stack the list into a time series of average correlation across sectors over time.

# the following computes average of the upper traingle correlation matrix elements
COR.mean <- sapply(COR.list, function(x) mean(x[upper.tri(x)])  )
summary(COR.mean)
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max.    NA's 
0.02625 0.45288 0.60169 0.57255 0.70916 0.94313      24 

To retrieve back into a time series object, following trick should serve well:

library(lubridate)
names(COR.mean) <- date(COR)
COR.mean <- as.xts(COR.mean) 
plot(COR.mean)

Note that, in order to transform a numerical vector into a time series, I label the values with the corresponding date and, then, set it as an xts object, whereas the lubridate is an extremely useful package to handle date formats.

LS0tCnRpdGxlOiAiVGlwIG9mIHRoZSBNb250aDogQ29ycmVsYXRpb24gT3ZlciBUaW1lIgojb3V0cHV0OiBybWFya2Rvd246OmdpdGh1Yl9kb2N1bWVudApvdXRwdXQ6CiAgaHRtbF9ub3RlYm9vazogZGVmYXVsdAogIHBkZl9kb2N1bWVudDogZGVmYXVsdAphdXRob3I6IE1hamVlZCBTaW1hYW4KZGF0ZTogTWF5IDE4LCAyMDE4CmZpZ193aWR0aDogNTAKLS0tCkluIG15IGVhcmxpZXIgcG9zdCBmcm9tIE1hcmNoIDIwMTgsIEkgZGlkIGludHJvZHVjZSB0aGUgYHJvbGxhcHBseWAgZnVuY3Rpb24gdGhhdCBleGVjdXRlcyBhIGZ1bmN0aW9uIG9uIGEgcm9sbGluZyB3aW5kb3cgYmFzaXMuIFdoaWxlIHRoaXMgZnVuY3Rpb24gaXMgdmVyeSB1c2VmdWwsIGl0IG5lZWRzIGEgbGl0dGxlIG1vZGlmaWNhdGlvbiBmb3IgdXNlcnMgdG8gYXBwbHkgb3RoZXIgZ2VuZXJhbCBvcGVyYXRpb25zLgoKT3JpZ2luYWxseSwgSSBmYWNlZCB0aGlzIGlzc3VlIHdoZW4gSSB0cmllZCB0byBjb21wdXRlIHRoZSBjb3JyZWxhdGlvbiBtYXRyaXggYWNyb3NzIGRpZmZlcmVudCBhc3NldCByZXR1cm5zIG9uIGEgcm9sbGluZyB3aW5kb3cuIEZvciB0aGUgZGVtb25zdHJhdGlvbiwgbGV0J3MgY29uc2lkZXIgdGhlIHJldHVybnMgZm9yIGFsbCBzZWN0b3IgRVRGcyBleGNsdWRpbmcgcmVhbCBlc3RhdGU6CmBgYHtyIG1lc3NhZ2U9RkFMU0UsIHdhcm5pbmc9RkFMU0V9CmxpYnJhcnkocXVhbnRtb2QpCnYgPC0gYygiWExFIiwiWExVIiwiWExLIiwiWExCIiwiWExQIiwiWExZIiwiWExJIiwiWExWIiwiWExGIikKdDEgPC0gIjE5OTAtMDEtMDEiClAubGlzdCA8LSBsYXBwbHkodixmdW5jdGlvbih4KSBnZXQoZ2V0U3ltYm9scyh4LGZyb20gPSB0MSkpICkKUC5saXN0IDwtIGxhcHBseShQLmxpc3QsZnVuY3Rpb24oeCkgeFssNl0pClAgPC0gUmVkdWNlKG1lcmdlLFAubGlzdCkKbmFtZXMoUCkgPC0gdgpSIDwtIG5hLm9taXQoUC9sYWcoUCkgLSAxKQpgYGAKQnkgZGVmYXVsdCwgYHJvbGxhcHBseWAgZXhlY3V0ZXMgdGhlIGdpdmVuIGZ1bmN0aW9uIG9uIGVhY2ggdGltZSBzZXJpZXMgc2VwYXJhdGVseSBhbmQgcmV0dXJucyBhIHRpbWUgc2VyaWVzIG9iamVjdC4gRm9yIGluc3RhbmNlLCAKYGBge3IgbWVzc2FnZT1GQUxTRSwgd2FybmluZz1GQUxTRX0KdGFpbChyb2xsYXBwbHkoUiwyNSxtZWFuKSkKYGBgCnJldHVybnMgdGhlIDI1IG1vdmluZyBhdmVyYWdlIGZvciBlYWNoIG9uZSBzZXBhcmF0ZWx5LiBPbiB0aGUgb3RoZXIgaGFuZCwgaWYgSSB0cnkgdG8gY29tcHV0ZSB0aGUgbW92aW5nIGNvcnJlbGF0aW9uLCBpbnN0ZWFkLCBJIGdldCB0aGUgZm9sbG93aW5nCmBgYHtyIG1lc3NhZ2U9RkFMU0UsIHdhcm5pbmc9RkFMU0V9CnRhaWwocm9sbGFwcGx5KFIsMjUsY29yKSkKYGBgCndoaWNoIGNvbXB1dGVzIHRoZSBjb3JyZWxhdGlvbiB3aXRoIHRoZSBzYW1lIEVURiByYXRoZXIgdGhhbiBvdGhlciBFVEZzIC0gYXMgaXQgdHJlYXRzIGVhY2ggdGltZSBzZXJpZXMgc2VwYXJhdGVseS4gQXMgYSByZW1lZHksIG9uZSBzaG91bGQgYWRkIGBieS5jb2x1bW4gPSBGYCBhcmd1bWVudCB0byB0aGUgYHJvbGxhcHBseWAgZnVuY3Rpb24uIEluIHRoaXMgY2FzZSwgdGhlIGZ1bmN0aW9uIHJldHVybnMgYSB0aW1lIHNlcmllcyBgeHRzYCBvYmplY3QsIGhvd2V2ZXIsIHdpdGggJDkgXHRpbWVzIDkgPSA4MSQgY29sdW1ucywgd2hlcmUgZWFjaCBjb2x1bW4gY29ycmVzcG9uZHMgdG8gdGhlIHBhaXJ3aXNlIGNvcnJlbGF0aW9uIGJldHdlZW4gdGhlIDkgc2VjdG9yIEVURnMgcmF0aGVyIHRoYW4gYSBzcXVhcmVkIG1hdHJpeC4KYGBge3IgbWVzc2FnZT1GQUxTRSwgd2FybmluZz1GQUxTRX0KQ09SIDwtIHJvbGxhcHBseShSLDI1LGNvcixieS5jb2x1bW4gPSBGKQpkaW0oQ09SKQpjbGFzcyhDT1IpCmBgYApXaGF0IGxlZnQgdG8gYmUgZG9uZSBpcyB0byBzdGFjayB0aGVzZSB2ZWN0b3JzIGJhY2sgaW50byBhIGNvcnJlbGF0aW9uIG1hdHJpeCwgb25lIGZvciBlYWNoIHRpbWUgcGVyaW9kLiBUbyBkbyBzbywgSSB3aWxsIHJlZmVyIHRvIHRoZSBgcGx5cmAgcGFja2FnZS4gVGhlIGBwbHlyYCBwYWNrYWdlIGFsbG93cyB1c2VycyB0byB0YWtlIGFuIGFycmF5IChgYWApLCBhIGRhdGEgZnJhbWUgKGBkYCksIG9yIGEgbGlzdCAoYGxgKSwgZXhlY3V0ZSBhIGdpdmVuIGZ1bmN0aW9uIG92ZXIgdGhlIGdpdmVuIG9iamVjdCwgYW5kIG91dHB1dCB0aGUgcmVzdWx0cyBpbiBlaXRoZXIgZm9ybWF0LiBGb3Igb3VyIGNhc2UsIEkgd2lsbCBpbnB1dCB0aGUgdGltZSBzZXJpZXMgYENPUmAgb2JqZWN0IGFzIGFuIGFycmF5IGFuZCBvdXRwdXQgaXQgYXMgYSBsaXN0LCB3aGVyZSBlYWNoIGVsZW1lbnQgaW4gdGhlIGxpc3QgIGNvcnJlc3BvbmRzIHRvIHRoZSBtb3ZpbmcgY29ycmVsYXRpb24gbWF0cml4LiAKYGBge3IgbWVzc2FnZT1GQUxTRSwgd2FybmluZz1GQUxTRX0KbGlicmFyeShwbHlyKQpDT1IubGlzdCA8LSBhbHBseShDT1IsMSxmdW5jdGlvbih4KSBtYXRyaXgoeCxucm93ID0gbmNvbChSKSwgIGJ5cm93ID0gVCApKQpgYGAKVGhlIHNlY29uZCBhcmd1bWVudCBpbiB0aGUgYGFscGx5YCBzcGVjaWZpZXMgdGhlIG1hcmdpbiwgd2hlcmUgMSBpbmRpY2F0ZXMgdGhhdCB0aGUgZ2l2ZW4gZnVuY3Rpb24gdG8gYmUgZXhlY3V0ZWQgb3ZlciB0aGUgcm93cywgd2hpbGUgMiBzdGF0ZXMgdGhhdCBpdCBzaG91bGQgYmUgZXhlY3V0ZWQgb3ZlciB0aGUgY29sdW1ucyBpbnN0ZWFkLiBUaGUgdGhpcmQgYXJndW1lbnQsIHdoaWNoIHRha2VzIGEgZnVuY3Rpb24sIHN0YWNrcyBlYWNoIHJvdyBvZiB0aGUgYENPUmAgb2JqZWN0IGludG8gYSBzcXVhcmVkIG1hdHJpeC4gQXMgYSByZXN1bHQsIHdlIGhhdmU6CmBgYHtyIG1lc3NhZ2U9RkFMU0UsIHdhcm5pbmc9RkFMU0V9CnJvdW5kKENPUi5saXN0W1syNV1dLDIpCmBgYAp3aGljaCBpcyBpZGVudGljYWwgdG8gY29ycmVsYXRpb24gbWF0cml4IGNvbXB1dGVkIG92ZXIgdGhlIGZpcnN0IDI1IGRheXMgaW4gdGhlIGRhdGEKYGBge3IgbWVzc2FnZT1GQUxTRSwgd2FybmluZz1GQUxTRX0Kcm91bmQoY29yKFJbMToyNSxdKSwyKQpgYGAKRmluYWxseSwgb25lIGNhbiBlaXRoZXIga2VlcCB0aGUgcm9sbGluZyBjb3JyZWxhdGlvbiBtYXRyaXggaW4gYSBsaXN0IG9yIHRyYW5zZm9ybSBpdCBiYWNrIGEgdGltZSBzZXJpZXMgdXNpbmcgY2VydGFpbiBjb21wdXRhdGlvbiwgZS5nLiBjb25zdHJ1Y3QgcG9ydGZvbGlvIHdlaWdodHMgYW5kIGNvbXB1dGUgdGhlIG91dC1vZi1zYW1wbGUgcmV0dXJuIGFzIGEgdGltZSBzZXJpZXMuIEFzIGEgZmluYWxsIGRlbW9zdHJhdGlvbiwgSSB3aWxsIHNob3cgaG93IG9uZSBjYW4gc3RhY2sgdGhlIGxpc3QgaW50byBhIHRpbWUgc2VyaWVzIG9mIGF2ZXJhZ2UgY29ycmVsYXRpb24gYWNyb3NzIHNlY3RvcnMgb3ZlciB0aW1lLgpgYGB7ciBtZXNzYWdlPUZBTFNFLCB3YXJuaW5nPUZBTFNFfQojIHRoZSBmb2xsb3dpbmcgY29tcHV0ZXMgYXZlcmFnZSBvZiB0aGUgdXBwZXIgdHJhaW5nbGUgY29ycmVsYXRpb24gbWF0cml4IGVsZW1lbnRzCkNPUi5tZWFuIDwtIHNhcHBseShDT1IubGlzdCwgZnVuY3Rpb24oeCkgbWVhbih4W3VwcGVyLnRyaSh4KV0pICApCnN1bW1hcnkoQ09SLm1lYW4pCmBgYApUbyByZXRyaWV2ZSBiYWNrIGludG8gYSB0aW1lIHNlcmllcyBvYmplY3QsIGZvbGxvd2luZyB0cmljayBzaG91bGQgc2VydmUgd2VsbDoKYGBge3IgbWVzc2FnZT1GQUxTRSwgd2FybmluZz1GQUxTRSxmaWcuYWxpZ249ImNlbnRlciJ9CmxpYnJhcnkobHVicmlkYXRlKQpuYW1lcyhDT1IubWVhbikgPC0gZGF0ZShDT1IpCkNPUi5tZWFuIDwtIGFzLnh0cyhDT1IubWVhbikgCnBsb3QoQ09SLm1lYW4pCmBgYAoKTm90ZSB0aGF0LCBpbiBvcmRlciB0byB0cmFuc2Zvcm0gYSBudW1lcmljYWwgdmVjdG9yIGludG8gYSB0aW1lIHNlcmllcywgSSBsYWJlbCB0aGUgdmFsdWVzIHdpdGggdGhlIGNvcnJlc3BvbmRpbmcgZGF0ZSBhbmQsIHRoZW4sIHNldCBpdCBhcyBhbiBgeHRzYCBvYmplY3QsIHdoZXJlYXMgdGhlIGBsdWJyaWRhdGVgIGlzIGFuIGV4dHJlbWVseSB1c2VmdWwgcGFja2FnZSB0byBoYW5kbGUgZGF0ZSBmb3JtYXRzLiAKCgoKCg==