使用 numpy 改进 worldquant 101 alpha 因子的实现

Question

我试图实现 WorldQuant 发布的 101 个量化交易因子 ( https://arxiv.org/pdf/1601.00991.pdf )。

\n

一个典型的因素是处理股票的价格和数量信息以及时间维度和股票维度。以 alpha 因子 #4 为例：(-1 * Ts_Rank(rank(low), 9))。这是动量阿尔法信号。最低价是一组股票在一定时期内的最低价。Rank 是对每行 panel\xe2\x80\x99s 进行排名的横截面过程（时间快照）。Ts_Rank 是在指定窗口内移动_rank panel\xe2\x80\x99s 每列（一只股票）的时间序列过程。

\n

直观上，Pandas 数据框或NumPy矩阵应该适合 101 个 alpha 因子的实现。下面是我迄今为止使用NumPy得到的最佳实现。然而，性能太低了。在我的 Intel core i7 Windows 机器上，以 5000（交易日期）× 200（股票）矩阵作为输入运行 alpha #4 因子大约需要 45 秒。

\n

我还遇到了 DolphinDB，一个具有内置分析功能的时间序列数据库（https://www.dolphindb.com/downloads.html）。对于相同的因子 Alpha#4，DolphinDB 的运行时间仅为 0.04 秒，比 NumPy 快 1000倍版本快 1000 倍。然而，DolphinDB 是商业软件。有人知道更好的 python 实现吗？或者有什么技巧可以改进我当前的 python 代码以实现与 DolphinDB 相当的性能？

\n

Numpy实现（基于https://github.com/yli188/WorldQuant_alpha101_code）

\n

import numpy as np\n\ndef rankdata(a, method=\'average\', *, axis=None):\n    # this rankdata refer to scipy.stats.rankdata (https://github.com/scipy/scipy/blob/v1.9.1/scipy/stats/_stats_py.py#L9047-L9153)\n    if method not in (\'average\', \'min\', \'max\', \'dense\', \'ordinal\'):\n        raise ValueError(\'unknown method "{0}"\'.format(method))\n\n    if axis is not None:\n        a = np.asarray(a)\n        if a.size == 0:\n            np.core.multiarray.normalize_axis_index(axis, a.ndim)\n            dt = np.float64 if method == \'average\' else np.int_\n            return np.empty(a.shape, dtype=dt)\n        return np.apply_along_axis(rankdata, axis, a, method)\n\n    arr = np.ravel(np.asarray(a))\n    algo = \'mergesort\' if method == \'ordinal\' else \'quicksort\'\n    sorter = np.argsort(arr, kind=algo)\n\n    inv = np.empty(sorter.size, dtype=np.intp)\n    inv[sorter] = np.arange(sorter.size, dtype=np.intp)\n\n    if method == \'ordinal\':\n        return inv + 1\n\n    arr = arr[sorter]\n    obs = np.r_[True, arr[1:] != arr[:-1]]\n    dense = obs.cumsum()[inv]\n\n    if method == \'dense\':\n        return dense\n\n    # cumulative counts of each unique value\n    count = np.r_[np.nonzero(obs)[0], len(obs)]\n\n    if method == \'max\':\n        return count[dense]\n\n    if method == \'min\':\n        return count[dense - 1] + 1\n\n    # average method\n    return .5 * (count[dense] + count[dense - 1] + 1)\n\ndef rank(x):\n    return rankdata(x,method=\'min\',axis=1)/np.size(x, 1)\n\ndef rolling_rank(na):\n    return rankdata(na.transpose(),method=\'min\',axis=0)[-1].transpose()    \n\ndef ts_rank(x, window=10):\n    a_rolled = np.lib.stride_tricks.sliding_window_view(x, window,axis = 0)\n    return np.append(np.full([window-1,np.size(x, 1)],np.nan),rolling_rank(a_rolled),axis = 0)\n\n\ndef alpha004(data):\n    return -1 * ts_rank(rank(data), 9)\n\nimport time\n\n# The input is a 5000 by 200 matrix, where the row index represents trade date and the column index represents security ID. \ndata=np.random.random((5000, 200))\nstart_time = time.time()\nalpha004(data)\nprint("--- %s seconds ---" % (time.time() - start_time))\n\n--- 44.85099506378174 seconds ---\n\n\n

\n

DolphinDB 实施

\n

def WQAlpha4(low){\n    return -mrank(rowRank(low, percent=true), true, 9)\n}\n\n// The input is a 5000 by 200 matrix, where the row index represents trade date and the column index represents security ID.\nlow = rand(1000.0,5000:200);\ntimer WQAlpha4(low);\n\nTime elapsed: 44.036 ms (0.044s)\n\n

\n

Answer 1

这部分代码：

return np.apply_along_axis(rankdata, axis, a, method)

……会很慢。像这样的函数应用程序意味着更多的计算在 Python 中运行，而相对较少的计算在 C 中运行。

如果您愿意对排名函数的定义方式稍作更改，那么这里有一个更快的解决方案。具体来说，下面的代码相当于从 更改method='min'为method='ordinal'。在随机数测试数据集上，它在 95% 的情况下与您的方法一致，只有 1 % 的情况不同。

通过沿轴使用 argsort，numpy 可以完成整个计算，而无需使用 Python。

def rank(x):
    return (data.argsort(axis=1).argsort(axis=1) + 1) / np.size(x, 1)


def ts_rank(x, window=10):
    a_rolled = np.lib.stride_tricks.sliding_window_view(x, window, axis = 0)
    rolling_rank_fast = (a_rolled.argsort(axis=2).argsort(axis=2) + 1)[:, :, -1]
    # Fill initial window - 1 rows with nan
    initial_window = np.full([window-1,np.size(x, 1)],np.nan)
    return np.append(initial_window,rolling_rank_fast,axis = 0)


def alpha004(data):
    return -1 * ts_rank(rank(data), 9)

通过对此进行基准测试，我发现它的运行速度大约快了 100 倍。