title: EWMA Covariance Matrix in Pandas - Optimization
body: <p>I would like to calculate the EWMA Covariance Matrix from a DataFrame of stock price returns using Pandas and have followed the methodology in <a href="https://github.com/robertmartin8/PyPortfolioOpt" rel="nofollow noreferrer">PyPortfolioOpt</a>.</p>
  <p>I like the flexibility of using Pandas objects and functions but when the set of assets grows the function is becomes very slow:</p>
  <pre><code>import pandas as pd
  import numpy as np
  
  def ewma_cov_pairwise_pd(x, y, alpha=0.06):
      x = x.mask(y.isnull(), np.nan)
      y = y.mask(x.isnull(), np.nan)
      covariation = ((x - x.mean()) * (y - y.mean()).dropna()
      return covariation.ewm(alpha=0.06).mean().iloc[-1]
  
  def ewma_cov_pd(rets, alpha=0.06):
      assets = rets.columns
      n = len(assets)
      cov = np.zeros((n, n))
      for i in range(n):
          for j in range(i, n):
              cov[i, j] = cov[j, i] = ewma_cov_pairwise_pd(
                  rets.iloc[:, i], rets.iloc[:, j], alpha=alpha)
      return pd.DataFrame(cov, columns=assets, index=assets)
  </code></pre>
  <p>I would like to improve the speed of the code ideally while still using Pandas but the bottleneck is within the DataFrame.ewm() function which uses 90% of the calculation time.</p>
  <p><strong>If using this function was a binding constraint, what is the most efficient way of improving the speed at which the code runs?</strong>  I was considering taking a brute force approach and using concurrent.futures.ProcessPoolExecutor but perhaps there is a better solutions.</p>
  <pre><code>n = 100  # n is typically 2000
  rets = pd.DataFrame(np.random.normal(0, 1., size=(n, n)))
  cov_pd = ewma_cov_pd(rets)
  </code></pre>
  <p>The true time-series data can contain leading nulls and potentially missing values after that although the latter less likely.</p>
  <p><strong>Update I</strong></p>
  <p>A potential solution which leverages off the answer provided by Quang Hoang and produces the expected results in a far more reasonable time would be something similar to:</p>
  <pre><code>def ewma_cov_frame_qh(rets, alpha=0.06):
      weights = (1-alpha) ** np.arange(len(df))[::-1]
      normalized = (rets-rets.mean()).to_numpy()    
      out = (weights * normalized.T) @ normalized / weights.sum()
      return pd.DataFrame(out, index=rets.columns, columns=rets.columns)
  
  
  def ewma_cov_qh(rets, alpha=0.06):
      syms = rets.columns
      covar = pd.DataFrame(index=rets.columns, columns=rets.columns)
      delta = rets.isnull().sum(axis=1).shift(1) - rets.isnull().sum(axis=1)
      dates = delta.loc[delta != 0].index.tolist()
       
      for date in dates:
          frame = rets.loc[rets.index &gt;= date].dropna(axis=1, how='any')
          cov = ewma_cov_frame_qh(frame).reindex(index=syms, columns=syms)
          covar = covar.fillna(cov)
     
      return covar
  
  cov_qh = ewma_cov_qh(rets)
  </code></pre>
  <p>This violates the requirement that the underlying covariance is calculated using the native Pandas/Numpy functions and calculation time will depend on the number leading na's in the data set.</p>
  <p><strong>Update II</strong></p>
  <p>A potential improvement on the above which uses (a naive implementation of) multiprocessing and improves the calculation time by a further 42.5% on my machine is listed below:</p>
  <pre><code>from concurrent.futures import ProcessPoolExecutor, as_completed
  from functools import partial
      
  def ewma_cov_mp_worker(date, rets, alpha=0.06):
      syms = rets.columns
      frame = rets.loc[rets.index &gt;= date].dropna(axis=1, how='any')
      return ewma_cov_frame_qh(frame, alpha=alpha).reindex(index=syms, columns=syms)
  
  
  def ewma_cov_mp(rets, alpha=0.06):
      covar = pd.DataFrame(index=rets.columns, columns=rets.columns)
      delta = rets.isnull().sum(axis=1).shift(1) - rets.isnull().sum(axis=1)
      dates = delta.loc[delta != 0].index.tolist()
  
      func = partial(ewma_cov_mp_worker, rets=rets, alpha=alpha)
      covs = {}
  
      with ProcessPoolExecutor(max_workers=6) as exec:
          future_to_date = {exec.submit(func, date): date for date in dates}
          covs = {future_to_date[future]: future.result() for future in as_completed(future_to_date)}
  
      for date in dates:
          covar.fillna(covs[date], inplace=True)
  
      return covar
  </code></pre>
  <p>[I have not added as answer as not addressed the original question and I am optimistic there is a better solution.]</p>