Canonical sectors and evolution of firms in
the US stock markets

Company names, tickers, listed-sectors and market caps of US-based firms used in this analysis were obtained from Scottrade^® (Scottrade^®, 2015). Daily closing prices adjusted for stock splits and dividend issues were obtained from Yahoo^® Finance (Yahoo!^® Finance, 2015). The rare cases of missing prices in the time series were replaced with linearly interpolated values. A brief summary of listed sectors and number of companies in each is provided in Table 1 and a full list of company names, tickers, market caps and listed-sector info is available on the companion website (Chachra et al., 2013).

A variety of factorization algorithms have been developed in recent years for dimensional reduction, classification or clustering. Examples include archetypal analysis (AA) (Cutler and Breiman, 1994), heteroscedastic matrix factorization (Tsalmantza and Hogg, 2012), binary matrix factorization (Zhang et al., 2007), K-means clustering (Ding and He, 2004), simplex volume maximization (Thurau et al., 2010), independent component analysis (Hyvärinen and Oja, 2000), non-negative matrix factorization (NMF) (Lee and Seung, 1999; Wang and Zhang, 2013) and its variants such as the semi- and convex-NMF (Ding et al., 2010), convex hull NMF (Thurau et al., 2011) and hierarchical convex NMF (Kersting et al., 2010), among others. Each method has a unique interpretation (Li and Ding, 2006) and therefore, a successful application of any of these methods is contingent upon the underlying structure of the data.

The hyper-tetrahedral structure of log price returns seen in our analysis motivates a decomposition so that each stock’s return is a weighted mixture of canonical sectors, constrained to lie in the convex hull of the data. Hence we employ AA factorization which is defined as:

Columns of $R_{ts}C_{sf}=E_{tf}$ are the emergent sector time series (basis vectors) representing the $n$ corners of the hyper-tetrahedron, and $W_{fs}$ are the participation weights ($W_{fs}\geq 0$) in sector $f$ so that $\sum_{f}W_{fs}=1$ for each stock $s$. The sector matrix $E_{tf}$ is within the convex hull ($C>0$, $\sum_{s}C_{sf}=1$) of the data $R_{ts}$. It can be found by either minimizing the squared error with convex constraints in factorization as originally proposed (Cutler and Breiman, 1994), or by making a convex hull of the dataset and choosing one or more of its vertices to be basis vectors, or by making a convex hull in low-dimensions and choosing one or more of its vertices to be basis vectors (Thurau et al., 2009), or by minimizing after initializing with candidate archetypes that are guaranteed to lie in the minimal convex set of the data (Mörup and Hansen, 2012). The columns of the $C$ matrix are shown in Fig. 6.

Numerical computations were performed using an in-house Python language implementation of the principal convex hull analysis (PCHA) algorithm as described in (Mörup and Hansen, 2012). For the full dataset, the factorization $R=EW$, with $E=RC$ as defined in Eqn. 1 converged in 35 iterations to a predefined tolerance value of $\Delta_{SSE}<10^{-7}$, where $\Delta_{SSE}$ is the average difference in the sum of squared error per matrix element in $R-EW$ from one iteration to the next. The resulting columns of $E_{tf}$ are shown in Fig. 7 (top row). Annualized cumulative log returns are obtained by summing rows of $E_{tf}$:

The time series $Q_{f}(\tau)$ are shown in Fig. 3 and the middle row of Fig. 7. Weights $W_{fs}$ for selected stocks are shown in Fig. 2, the remainder are available on the companion website (Chachra et al., 2013). In each canonical sector $f$, the component of weights for companies are shown in Fig. 8.

The analysis of evolving sector weights was performed similarly, but with a sliding Gaussian time window. We decomposed the local normalized log returns for each stock into the canonical sectors determined from the entire time series. Each column (time series) of the returns matrix $R_{ts}$ was multiplied with a Gaussian, $G_{\mu}(\tau)=\exp(-(\tau-\mu)^{2}/(2\times 250^{2}))$ of standard deviation $250$ centered at $\mu$ to obtain $R^{\mu}_{ts}$. We use $C_{s^{\prime}f}$ found using the full dataset (Eqn. 1) (corresponding to keeping the sector-defining simplex corners fixed). $R^{\mu}_{ts}$ is factorized to obtain new weights $W^{\mu}_{fs}$ that describe sector decomposition of stocks in that period focused at $t=\mu$: $R^{\mu}=R^{\mu}_{ts^{\prime}}C_{s^{\prime}f}W^{\mu}_{fs}$. $\mu$ is increased in steps of 50 starting at $\mu=0$ and ending at $\mu=5000$, and $W^{\mu}$ is calculated at each $\mu$ with the corresponding $R^{\mu}$. These results are plotted in Fig. 5 for a select group of companies; the remainder are available on the companion website (Chachra et al., 2013).

To address the challenge of distinguishing signal from noise in the evolving sector weights, we emulated the effect of noise for each of the companies from Fig. 5. For each of these companies, we took its sector weights, $\vec{\omega_{f}}$, and multiplied by $E_{tf}$ to obtain a time series for the company with weights that are constant in time. We then added gaussian random noise with standard deviation one and replaced these companies by this simulated data. Fig. 9 shows the comparison between the real flows and the simulated constant data with noise added. General features are shown to be signal while small fluctuations are consistent with noise.

It is often the case with large datasets that the effective dimensionality of the data space is much lower when one filters out the noise. Of the many dimensional reduction methods, the most commonly used is singular value decomposition (SVD) (Press et al., 2007), a deterministic matrix factorization. We discuss SVD in more detail in order to draw a contrast with previous SVD results, and to apply it for quantifying the explainable variation in the returns data.

An SVD of $R_{ts}$ is a matrix factorization (Press et al., 2007) $R_{ts}=U_{tf}\Sigma_{ff^{\prime}}V^{T}_{f^{\prime}s}$ such that matrices $U$ and $V$ are orthogonal; $\Sigma$ is a diagonal matrix of “singular values”. If the goal were purely rank-reduction, $n$ entries of $\Sigma$ chosen to lie above “noise threshold” are retained and the rest truncated so that $0\leq f,\ f^{\prime}\leq n$. This effectively reduces the dimension of $R$ to $n$. The choice of $n$ can be informed by the distribution of singular values as discussed later. The rows of $V^{T}$ are precisely the eigenvectors of the stock-stock returns correlation matrix, $\xi_{ss^{\prime}}\sim R_{st}^{T}R_{ts}$. It was previously reported that some components of the stiff eigenvectors of this stock-stock correlation matrix loosely corresponded to firms belonging to the same conventionally identified business sector (Plerou et al., 2002) (but see Fig. 10).

After normalizing the log returns, the returns matrix $R$ has entries of unit variance. If the entries were uncorrelated random variables drawn from a standard normal distribution, their singular values (which are also the positive square roots of the eigenvalues of $R^{T}R$) would be described by Wishart statistics (Mehta, 2004). The Wishart ensemble for a matrix of size $\alpha\times\beta$ predicts a distribution of singular values with a characteristic shape (Mehta, 2004), bounded for large matrices by $\sqrt{\alpha}\pm\sqrt{\beta}$. Comparing the stock correlations with Wishart statistics has been previously used to filter noise from financial datasets (Laloux et al., 1999). As shown in Fig. 11, most singular values of the returns matrix $R$ lie in the bulk below the bound set by the Wishart ensemble, whereas only $\sim$20 fall outside that cutoff (The singular value bounds of a random Gaussian rectangular matrix of size $\alpha\times\beta$ can be shown to be $\sqrt{\alpha}\pm\sqrt{\beta}$ for large matrices.) Historically, this has served as indication that singular values within the bulk correspond to noise (Laloux et al., 1999). Recently, however, much progress has been made in the development of techniques to extract signal from the bulk (Burda et al., 2004, 2006; Livan et al., 2011). Our method does not claim to capture this information. Rather, we measure its ability to capture variation in the data above the cutoff by means of random matrix theory explainable variation as defined in section F. The largest singular value of $R_{ts}$ corresponds to what we will refer to as the “market mode” as this represents overall simultaneous rise and fall of stocks. In the analysis presented in this paper, this mode has been filtered from the returns matrix by projecting the $R$ matrix into the subspace spanned by all non-market mode eigenvectors. This is nearly equivalent to filtering the market mode using simple linear regression (as done commonly (Plerou et al., 2002)), although more convenient.

The emergent low-dimensional, hyper-tetrahedral (simplex) structure of stock price returns can be seen by projecting the dataset into stiff “eigenplanes”. Eigenplanes are formed by pairs of right singular vectors from a SVD. Here, we construct an SVD of the simplex corners, $E_{tf}=X_{tk}YZ^{T}_{kf}$; simplex corners are mapped to columns of $YZ^{T}$ because $YZ_{kf}^{T}=X_{kt}^{T}E_{tf}$ (in other words, $X^{T}_{kt}$ is a projection operator). The plots in Fig. 12 are the projections of the dataset, $X^{T}_{kt}R_{ts}=v_{ks}$. The rows of $v$ taken in pairs form the axes of the projections in Fig. 1 and Fig. 12. With those plots, it becomes clear that the eigenplanes represent projections of a simplex-like data into two-dimensions. Secondly, we note that the simplex structure becomes less clear as one looks at planes corresponding to smaller singular value directions; the signal eventually becomes buried in the noise.

Similarly, the results of the factorization can be seen in eigenplanes from the SVD of $E_{tf}W_{sf}=L_{tk}MN^{T}_{ks}$. These results (rows of $MN^{T}_{ks}$) are shown in Fig. 13, where we notice that the data is now perfectly resides in simplex region as expected due to constraints.

We measured the goodness of the returns decomposition $R=EW$ by measuring the coefficient of determination ($r^{2}$) as follows:

Here, SSE denotes the sum of square errors $||R-EW||_{F}^{2}$, and SST is the total sum of squares $||R||^{2}_{F}$. This is also known as the proportion of variance explained (PVE). For the factorization of the full dataset, normalized with the market mode removed, the calculated $r^{2}$ value is $11.1\%$. The SVD of $R$ with singular values shown in Fig. 11 provides a convenient way to put this number in context for the returns dataset. Only 20 singular values (excluding the market mode) were above the cut-off that was predicted by random matrix theory for a matrix of purely random Gaussian entries. For any matrix $M$ with elements $m_{ij}$, the norm $||M||^{2}_{F}=\sum_{i,j}m_{ij}^{2}=\sum_{i}s_{i}^{2}$, where $s_{i}$ are the singular values (Press et al., 2007). Thus, the fraction of intrinsic variation in $R$ above the cutoff is the sum of squares of the 20 singular values (not including market mode) divided by SST, $\sum_{i=1}^{i=20}s_{i}^{2}/||R||^{2}_{F}=19.8\%$. Therefore, as a first approximation, the factorization explains $11.1/19.8=56\%$ of the random matrix theory (RMT) explainable variation.

For reference we provide the RMT explainable variation for the factor decomposition of Fama and French, the classification by Scottrade^®, and the top 8 singular vectors given by SVD. The percentage of the RMT explainable variation for different numbers of factors compared to the 3 factor decomposition of Fama and French is shown in Table 2. Fama and French have the benefit of allowing factors to have positive or negative weights. In order to compare with another non-negative decomposition, we fix the weight matrix according to the Scottrade^® labels and run archetypal analysis for this $n=14$ factor version. The $r^{2}$ value for this decomposition is $10.7\%$ with a corresponding RMT explainable variance of $54.2\%$ compared to $56\%$ for our 8 factors. For completeness, we also note that if $R$ is rank-reduced to the eight stiffest components found by SVD (not including market mode), then the factorization explains $85\%$ of the the RMT explainable variation in $R$ with overall results in good accord with the analysis presented here. This implies that sector decomposition information was already contained in the stiff modes from the SVD of $R$, however SVD is not the appropriate tool for the decomposition. Fig. 14 further shows that our unsupervised 3-factor decomposition appears quite distinct from Fama and French’s hand-created one.

It is an open problem to determine the effective dimensionality (optimal rank) of a general dataset (matrix). One could select among models of different dimensions using statistical tests such as the $r^{2}$ discussed above, or information theory based criteria such as Akaike Information Criterion (AIC) or the Bayesian Information Criterion (BIC), but the choice of the selection criterion is itself generally made on an ad hoc basis. Therefore, a direct observation of the comprehensibility of results is often the most reliable criterion. In the dataset used for analysis described here, a factorization with $n>8$ yielded results where both the emergent time series $E_{tf}$ and weights in $W_{fs}$ showed qualitative signs of overfitting. For example, with $n=9$ the results were in good agreement with $n=8$ except for an additional resulting sector involving participation from only 11 seemingly unrelated stocks (Table 3 and Fig. 4). The high-level results of factorization with different values of $n$ may be explored in a number of ways, several of which are described below.

The matrix $C_{sf}$ in decomposition $R=RCW$ represents how returns $R$ of stocks $s$ must be combined to make canonical sector returns $E_{tf}=R_{ts}C_{sf}$. Since a canonical sector is defined as a combination of stocks, an investment in the sector $f$ can made via buying a basket of constituent stocks $s$ in proportions given by $C_{sf}$ or through an index $I_{tf}$:

where, $p$ are stocks prices suitably weighted by market cap or other divisor as common practice for common indices (Tagiliani, 2009). An unweighted index of this kind is shown in the bottom row of Fig. 7 for results corresponding to the analysis described in this paper. Conversely, a pre-defined basket of stocks such as the S&P 500^® can be unbundled to find its exposure to the canonical sectors. With an investment strategy employing longs and shorts at the same time in correct proportions, it is conceivable to invest in, for example, the c-tech component of S&P 500^®.

The desirable features of an index include completeness, objectivity and investability (Pastor et al., 2013). The c-indices constructed using the ideas outlined here would not only be of value to investors through investment vehicles such as ETFs, Futures, etc., but also serve as important economic indicators.