このページの内容は最新ではありません。最新版の英語を参照するには、ここをクリックします。

Tikhonov 正則化を使用した実数最小二乗行列解の固定小数点型の判別

この例では、関数fixed.realQRMatrixSolveFixedpointTypesを使用して、次に示す実数最小二乗行列方程式の解の固定小数点型を解析的に判別する方法を説明します。

$$\left[\begin{array}{c}\lambda I_n\\ A\end{array}\right]X=\left[\begin{array}{c}{0}_{n,p}\\ B\end{array}\right],$$

ここで、$$A$$ は $$m$$ 行 $$n$$ 列の行列 ($$m\ge n$$)、$$B$$ は $$m$$ 行 $$p$$ 列、$$X$$ は $$n$$ 行 $$p$$ 列、$$I_n=\text{eye}(n)$$、$$0_{n,p}=\text{zeros}(n,p)$$、$$\lambda$$ は正則化パラメーターです。

最小二乗解は次のとおりです。

$$X_{LS}=({\lambda }^2{I}_n+A^{T}A{)}^{-1}{A}^{T}B$$

ただし、二乗や逆数を使用せずに計算されます。

システム パラメーターの定義

この例の行列属性とシステム パラメーターを定義します。

m は、行列 A および B の行数です。ビームフォーミングや方向探知などの問題では、m は積分するサンプル数に対応します。

m = 300;

n は、行列 A の列数と行列 X の行数です。最小二乗問題では、m は n より大きく、通常、m は n より大幅に大きくなります。ビームフォーミングや方向探知などの問題では、n はセンサー数に対応します。

n = 10;

p は、行列 B および X の列数です。p 個の右辺をもつ系の同時求解に対応します。

p = 1;

この例では、行列 A のランクを列数未満に設定します。ビームフォーミングや方向探知などの問題では、$$\text{rank}(A)$$ はセンサー アレイに作用する信号の数に対応します。

rankA = 3;

precisionBits は、行列解に必要な精度のビット数を定義します。この値はシステム要件に従って設定します。

precisionBits = 32;

正則化パラメーターを小さい正の値にすると、問題の調整を向上させ、推定値の分散を小さくすることができます。バイアスがありますが、推定値の分散が小さいと、多くの場合、最小二乗推定値と比べて平均二乗誤差が小さくなります。

regularizationParameter = 0.01;

この例では、実数値の行列 A および B は、それらの要素の大きさが 1 以下となるように構成されています。これらの値は各自のシステム要件によって定義されます。これらの内容が不明な場合に、A と B がシステムへの固定小数点入力である場合は、関数upperboundを使用して A と B の固定小数点型の上限を判別できます。

max_abs_A は、A の最大の大きさ要素に対する上限です。

max_abs_A = 1;

max_abs_B は、B の最大の大きさ要素に対する上限です。

max_abs_B = 1;

熱ノイズの標準偏差は、システム パラメーターである熱ノイズ パワーの平方根です。適切に設計されたシステムの量子化レベルは熱ノイズより低くなります。ここでは、thermalNoiseStandardDeviation を $$-50$$ dB のノイズ パワーと等価に設定します。

thermalNoiseStandardDeviation = sqrt(10^(-50/10))

thermalNoiseStandardDeviation = 0.0032

量子化ノイズの標準偏差は、必要な精度のビット数の関数です。これを計算するには、fixed.realQuantizationNoiseStandardDeviationを使用します。これが thermalNoiseStandardDeviation 未満であることを確認します。

quantizationNoiseStandardDeviation = fixed.realQuantizationNoiseStandardDeviation(precisionBits)

quantizationNoiseStandardDeviation = 6.7212e-11

固定小数点型の計算

この例では、設計されたシステム行列 $$A$$ がフル ランクでなく (対象の信号の数が行列 $$A$$ の列数より少ない)、測定されたシステム行列 $$A$$ の加法性熱ノイズが量子化ノイズより大きいと仮定します。この加法性ノイズにより、測定された行列 $$A$$ はフル ランクになります。

$${\sigma }_{\text{noise}}={\sigma }_{\text{thermal}\text{noise}}$$ を設定します。

noiseStandardDeviation = thermalNoiseStandardDeviation;

fixed.realQRMatrixSolveFixedpointTypesを使用して固定小数点型を計算します。

T = fixed.realQRMatrixSolveFixedpointTypes(m,n,max_abs_A,max_abs_B,...
    precisionBits,noiseStandardDeviation,[],regularizationParameter)

T = struct with fields:
    A: [0x0 embedded.fi]
    B: [0x0 embedded.fi]
    X: [0x0 embedded.fi]

T.A は、$$\left[\begin{array}{c}\lambda I_n\\ A\end{array}\right]$$ をオーバーフローしないようにインプレースで $$R=Q^T\left[\begin{array}{c}\lambda I_n\\ A\end{array}\right]$$ に変換するために計算された型です。

T.A

ans = 

[]

          DataTypeMode: Fixed-point: binary point scaling
            Signedness: Signed
            WordLength: 39
        FractionLength: 32

T.B は、$$\left[\begin{array}{c}{0}_{n,p}\\ B\end{array}\right]$$ をオーバーフローしないようにインプレースで $$C=Q^T\left[\begin{array}{c}{0}_{n,p}\\ B\end{array}\right]$$ に変換するために計算された型です。

T.B

ans = 

[]

          DataTypeMode: Fixed-point: binary point scaling
            Signedness: Signed
            WordLength: 39
        FractionLength: 32

T.X は、オーバーフローする可能性が低くなるような解 $$X=\left[\begin{array}{c}\lambda I_n\\ A\end{array}\right]\backslash \left[\begin{array}{c}{0}_{n,p}\\ B\end{array}\right]$$ を求めるために計算された型です。

T.X

ans = 

[]

          DataTypeMode: Fixed-point: binary point scaling
            Signedness: Signed
            WordLength: 44
        FractionLength: 32

指定した型を使用して行列方程式を解く

B が A の範囲内にあり、rankA=rank(A) となるような乱数行列 A および B を作成します。ランダム測定ノイズを A に追加してこれをフル ランクにしますが、この影響は解にも及び、B は A の範囲に近いだけになります。

rng('default');
[A,B] = fixed.example.realRandomLeastSquaresMatrices(m,n,p,rankA);
A = A + fixed.example.realNormalRandomArray(0,noiseStandardDeviation,m,n);

入力をfixed.realQRMatrixSolveFixedpointTypesで判別された型にキャストします。固定小数点に量子化することは、ランダム ノイズ [4,5] を追加することと等価です。

A = cast(A,'like',T.A);
B = cast(B,'like',T.B);

fiaccelを使用して MATLAB 実行可能 (MEX) 関数を生成することで、関数fixed.qrMatrixSolveを高速化します。

fiaccel fixed.qrMatrixSolve -args {A,B,T.X,regularizationParameter} -o qrRealMatrixSolve_mex

出力型 T.X を指定し、QR 法を使用して固定小数点 $$X=A\backslash B$$ を計算します。

X = qrRealMatrixSolve_mex(A,B,T.X,regularizationParameter);

出力の精度の検証

既定の倍精度浮動小数点値を使用した固定小数点出力と MATLAB からの出力の間の相対誤差が小さいことを確認します。

$$X_{\text{double}}=\left[\begin{array}{c}\lambda I_n\\ A\end{array}\right]\backslash \left[\begin{array}{c}{0}_{n,p}\\ B\end{array}\right]$$

A_lambda = double([regularizationParameter*eye(n);A]);
B_0 = [zeros(n,p);double(B)];
X_double = A_lambda\B_0;
relativeError = norm(X_double - double(X))/norm(X_double)

relativeError = 4.4740e-06

このファイルでは mlint の警告は非表示にします。

%#ok<*NASGU>
%#ok<*ASGLU>

参考

fixed.realQRMatrixSolveFixedpointTypes | Real Burst Matrix Solve Using QR Decomposition | Real Partial-Systolic Matrix Solve Using QR Decomposition