インテル® MKL クックブック

金融市場のデータストリームにおけるノイズ・フィルタリング

目的

一部の株式の価格変動が大規模な株式ポートフォリオのほかの株式の価格変動にどのように影響を及ぼすかを特定する。

ソリューション

データ全体の依存関係を表す相関行列を、信号行列とノイズ行列の 2 つの成分に分割します。信号行列から、株式間の依存関係を正確に評価することができます。アルゴリズム ( [Zhang12]、[Kargupta02]) は、累積データの相関行列の固有値を調べることにより、有用な情報からノイズ成分を分割する固有状態ベースのアプローチを取ります。

インテル® MKL サマリー統計は、ストリーミング・データの相関行列を計算する機能を提供します。インテル® MKL LAPACK には、さまざまな特性や格納形式の対称行列の固有値と固有ベクトルを計算する計算ルーチン群が含まれています。

オンライン・ノイズ・フィルタリング・アルゴリズムの手順を次に示します。

  1. λmin と λmax (ノイズ固有状態の間隔の境界) を計算します。

  2. データストリームから新しいデータブロックを取得します。

  3. 最新のデータブロックを使用して相関行列を更新します。

  4. 間隔 [λmin, λmax] に属する相関行列の固有値を検索して、ノイズ成分を定義する固有値と固有ベクトルを計算します。

  5. ステップ 4 で計算した固有値と固有ベクトルを組み合わせてノイズ成分の相関行列を計算します。

  6. 相関行列全体からノイズ成分を取り除いて信号成分の相関行列を計算します。さらにデータがある場合は、ステップ 2 に戻ります。

ソースコード: サンプル (http://software.intel.com/en-us/mkl_cookbook_samples) の nf フォルダーを参照してください。

初期化

相関解析タスクとパラメーターを初期化します。 

VSLSSTaskPtr task;
double *x, *mean, *cor;
double W[2];
MKL_INT x_storage, cor_storage;

...

scanf("%d", &m);            // ブロックの観測数
scanf("%d", &n);            // 株式の数 (タスクの次元)

...

/* メモリーを割り当て */
nfAllocate(m, n, &x, &mean, &cor, ...);


/* サマリー統計タスク構造を初期化 */
nfInitSSTask(&m, &n, &task, x, &x_storage, mean, cor, &cor_storage, W);
...

/* メモリーの割り当て */
void nfAllocate(MKL_INT m, MKL_INT n, double **x, double **mean, double **cor,
               ...)
{
    *x = (double *)mkl_malloc(m*n*sizeof(double), ALIGN);
    CheckMalloc(*x);

    *mean = (double *)mkl_malloc(n*sizeof(double), ALIGN);
    CheckMalloc(*mean);

    *cor = (double *)mkl_malloc(n*n*sizeof(double), ALIGN);
    CheckMalloc(*cor);
...
}
/* サマリー統計タスク構造を初期化 */
void nfInitSSTask(MKL_INT *m, MKL_INT *n, VSLSSTaskPtr *task, double *x,
                  MKL_INT *x_storage, double *mean, double *cor,
                  MKL_INT *cor_storage, double *W)
{
    int status;

    /* VSL サマリー統計タスクを作成 */
    *x_storage = VSL_SS_MATRIX_STORAGE_COLS;
    status = vsldSSNewTask(task, n, m, x_storage, x, 0, 0);
    CheckSSError(status);

    /* タスクの重みのレジスター配列 */
    W[0] = 0.0;
    W[1] = 0.0;
    status = vsldSSEditTask(*task, VSL_SS_ED_ACCUM_WEIGHT, W);
    CheckSSError(status);

    /* 相関行列計算にフルストレージを使用する
       タスク・パラメーターの初期化 */
    *cor_storage = VSL_SS_MATRIX_STORAGE_FULL;
    status = vsldSSEditCovCor(*task, mean, 0, 0, cor, cor_storage);
    CheckSSError(status);
}

計算

データの各ブロックについてノイズ・フィルタリングのステップを実行します。 

/* ノイズ成分を定義するしきい値を設定 */
sqrt_n_m = sqrt((double)n / (double)m);
lambda_min = (1.0 - sqrt_n_m) * (1.0 - sqrt_n_m);
lambda_max = (1.0 + sqrt_n_m) * (1.0 + sqrt_n_m);

...
/* データブロックをループ */
for (i = 0; i < n_block; i++)
{
    /* データの次の部分を読む */
    nfReadDataBlock(m, n, x, fh);

    /* "信号" と "ノイズ" 共分散の評価を更新 */
    nfKernel(m, n, lambda_min, lambda_max, x, cor, cor_copy,
             task, eval, evect, work, iwork, isuppz,
             cov_signal, cov_noise);
}
...

void nfKernel(…)
{
...

    /* FAST メソッドを使用して相関行列の評価を更新 */
    errcode = vsldSSCompute(task, VSL_SS_COR, VSL_SS_METHOD_FAST);
    CheckSSError(errcode);

...
    /* 相関行列の固有値と固有ベクトルを計算 */
    dsyevr(&jobz, &range, &uplo, &n, cor_copy, &n, &lmin, &lmax,
           &imin, &imax, &abstol, &n_noise, eval, evect, &n, isuppz,
           work, &lwork, iwork, &liwork, &info);

    /* 共分散行列の "信号" と "ノイズ" 部分を計算 */
    nfCalculateSignalNoiseCov(n, n_signal, n_noise,
        eval, evect, cor, cov_signal, cov_noise);
}

...
static int nfCalculateSignalNoiseCov(int n, int n_signal, int n_noise,
        double *eval, double *evect, double *cor, double *cov_signal,
        double *cov_noise)
{
    int i, j, nn;

    /* SYRK パラメーター */
    char uplo, trans;
    double alpha, beta;

    /* 共分散行列の "ノイズ" 部分を計算 */
    for (j = 0; j < n_noise; j++) eval[j] = sqrt(eval[j]);

    for (i = 0; i < n_noise; i++)
        for (j = 0; j < n; j++)
            evect[i*n + j] *= lambda[i];

    uplo  = 'U';
    trans = 'N';
    alpha = 1.0;
    beta  = 0.0;
    nn = n;

    if (n_noise > 0)
    {
        dsyrk(&uplo, &trans, &nn, &n_noise,  &alpha, evect, &nn,
              &beta, cov_noise, &nn);
    }
    else
    {
        for (i = 0; i < n*n; i++) cov_noise[i] = 0.0;
    }

    /* 共分散行列の "信号" 部分を計算 */
    if (n_signal > 0)
    {
        for (i = 0; i < n; i++)
            for (j = 0; j <= i; j++)
                cov_signal[i*n + j] = cor[i*n + j] - cov_noise[i*n + j];
    }
    else
    {
        for (i = 0; i < n*n; i++) cov_signal[i] = 0.0;
    }

    return 0;
}

リソースの解放

タスクを削除し、関連付けられたリソースを解放します。 

errcode = vslSSDeleteTask(task);
CheckSSError(errcode);
MKL_Free_Buffers();

使用するルーチン

タスク

ルーチン

説明

サマリー統計タスクを初期化して解析用のオブジェクト (データセット、サイズ (変数の数と観測の数)、格納形式)) を定義する

vsldSSNewTask

新しいサマリー統計タスク・ディスクリプターを作成して初期化します。

相関行列を保持するメモリーを指定する

vsldSSEditCovCor

共分散/相関/クロス積パラメーターのポインターを変更します。

処理した観測の重みの合計を保持する 2 要素の配列を指定する (データストリームの評価の正確な計算に必要)

vsldSSEditTask

タスク・ディスクリプターの入出力パラメーターのアドレスを変更します。

計算タイプ (VSL_SS_COR) と計算メソッド (VSL_SS_METHOD_FAST) を指定して主計算ドライバーを呼び出す

vsldSSCompute

サマリー統計評価を計算します。

タスクに関連付けられているリソースの割り当てを解除する

vslSSDeleteTask

タスク・オブジェクトを破棄してメモリーを解放します。

相関行列の固有値と固有ベクトルを計算する

dsyevr

選択した固有値を (オプションで、Relatively Robust Representation アルゴリズムを使用して実対称行列の固有値を) 計算します。

対称ランク-k 更新する

dsyrk

対称ランク-k 更新を行います。

説明

アルゴリズムのステップ 4 では、対称行列の固有値問題を解きます。オンライン・ノイズ・フィルタリング・アルゴリズムを使用するには、データのノイズを定義する、事前定義間隔 [λmin, λmax] に属する固有値を計算する必要があります。 LAPACK ドライバールーチン ?syevr は、この種の問題を解くデフォルトのルーチンです。 ?syevr インターフェイスで、呼び出し元は、固有値を検索する範囲の下限と表現のペア (このケースでは λmin と λmax) を指定できます。

固有ベクトルは直交行列 A を含む配列の列として返され、固有値は対角行列 Diag の要素を含む配列で返されます。 ノイズ成分の相関行列は、A*Diag*AT を計算して取得できます。 しかし、ここでは、2 つの一般的な行列乗算によってノイズ相関行列を構築する代わりに、1 つの対角行列乗算と 1 つのランク n 更新操作を使用してより効率良く計算しています。

ランク n 更新操作用に、インテル® MKL は BLAS 関数 ?syrk を提供しています。

親トピック: インテル® MKL クックブック