株式会社シグナルアナリシス技術資料

Subspaceの考え方を基にした妨害信号除去

株式会社シグナルアナリシス

関原謙介

本文書の著作権は株式会社シグナルアナリシスに帰属します．

妨害信号
磁気シールドルームの外で発生する妨害信号（環境のノイズ）
実験室の電子機器，建てや内のエレベータ，空調設備，外を走る車や電車等が発生する磁場．

磁気シールドルームの中で発生する妨害信号
生体起源の妨害磁場（生理学的ノイズ）：眼球運動や心臓の発生する磁場，（誘発脳活動が測定対象なら）自発的な脳活動による磁場（brain noise)など．
非生体起源の磁場：DBSやVNSなどの脳刺激装置の発生器からの磁場．

本文書の著作権は株式会社シグナルアナリシスに帰属します．

妨害信号除去のために知られている技術
環境のノイズ
次のスライドでまとめる
眼球運動や心臓の発生する磁場
ICA

SSP

DSSP（多分）

ブレインノイズ（臨床応用で問題となるか？）PrewhiteningPFA（ベイズ因子分析を発展させた手法）

脳刺激装置の発生器からの磁場
tSSS

DSSP

Adaptive beamforming
本文書の著作権は株式会社シグナルアナリシスに帰属します．

環境ノイズ除去法のまとめ別に取られたノイズデータを必要とする方法

ノイズの空間情報を利用

ノイズの時間情報を利用

名称

noise dataの取得法

Signal SpaceProjection(SSP)

empty-roomnoise data等を使用

Adaptive noise
cancelling (ANC)

リファレンスセン
サーデータ

Common subspace
projection (CSP)

リファレンスセン
サーデータ

ノイズデータを必要としない方法

ノイズの時間・空間情報を利用•••

ノイズ発生源の空間配置が時間的に変わらないことを仮定．
信号にも影響を与える可能性がある．

適切に配置されたリファレンスセンサーが必要．

利点・欠点・特徴など

名称
ノイズの空間情報を利用

利点・欠点

Signal Space
Separation (SSS)••••

ヘルメット型のMEGセンサーアレイにのみ有効．
高精度なセンサーキャリブレーションが必要．
エレクタから特許出願有り．
信号源に近接した妨害信号源では有効性低下．

Spatio-temporal
SSS (tSSS)••

ヘルメット型のMEGセンサーアレイにのみ有効．
エレクタから特許出願有り．

Dual signalsubspace
projection (DSSP)•••

任意のセンサーアレイに有効．
いくつかの数値パラメータを調節する必要有り．
SN比の比較的良い場合に有効性低下．
本文書の著作権は株式会社シグナルアナリシスに帰属します．

Signal space projection algorithm
妨害信号と関心信号の空間的特性の違いを利用
エレクタ-Neuromagが1990年代半ばに，Neuromag-MEG装置のmagnetometerセンサーのノイズを取るために用いた方法．

y(t )データモデル =yS (t ) + ε

計測データからどうやって信号とノイズを分離するか？
Q

信号ベクトル： yS (t ) = ∑ sq (t )lqq =1

は時刻t に応じていろいろな値を取るが，M 次元空間の中である決まった領域を占める．
この領域を信号部分空間と呼び，span {l1, , lQ } で表される．

l1, , lQ の線形和で表される全てのベクトルのなす集合の意味である

本文書の著作権は株式会社シグナルアナリシスに帰属します．

計測データからどうやって信号とノイズを分離するか？

y(t )データモデル：=ε  1ノイズベクトル： ε =   ε  M

yS (t ) + ε

ランダムにいろんな方向を向く

信号ベクトル：
ソース1個の場合

yS (t ) = s1(t )l1

このベクトルは l 方向をいつも向いている
1

ソース2個の場合=
yS (t ) s1(t )l1 + s2 (t )l2

このベクトルは l と l の線形和で表される
12

方向をいつも向いている

span {l1, l2 } と呼ぶ
本文書の著作権は株式会社シグナルアナリシスに帰属します．

ソースQ個の場合
Q

信号ベクトル： yS (t ) = ∑ sq (t )lq は span {l , , l1q =1Q} で表される領域を占める

この領域を信号部分空間とぶ．
信号部分空間が分かれば，その情報を利用してノイズと信号の分離ができそう？
１）どうやって信号部分空間を見つけるか？（信号源の位置と向きは未知であるのに）
２）信号部分空間を見つけられたとして，どうやってノイズと信号の分離を行うのか
リードフィールドベクトル

l1, , lQ は線形独立であるので，信号部分空間はQ 次元である．

本文書の著作権は株式会社シグナルアナリシスに帰属します．

どうやってノイズと信号の分離を行うのかもし，信号部分空間の（正規直交）基底ベクトル u1, , uQ が分かれば，
QyS (t ) = ∑ cq uq と表すことができる．
q =1

T

プロジェクター P = [ u1, , uQ ][ u1, , uQ ]をデータベクトルに乗じて，Py(t )= PyS (t ) + Pε= yS (t ) + Pε である．
TQ=
PyS (t ) [ u1, , uQ ][ u1, , uQ ] ∑=cq uq

Q=cq uq yS (t )∑=q 1=q 1

Pε

ノイズベクトルの信号部分空間内の成分で，非常に小さい
本文書の著作権は株式会社シグナルアナリシスに帰属します．

どうやって，信号部分空間の基底ベクトルu1, , uQ を求めるか？

もし，リードフィールドベクトル l1, , lQ が分かっているなら，l1, , lQ を直交化することにより u1, , uQ を求めることができる．
しかし， l , , l1

Qは未知量である．

l1, , lQ を知らずしていかにして u1, , uQ を求めるか？

本文書の著作権は株式会社シグナルアナリシスに帰属します．

l1, , lQ を知らずしていかにして u1, , uQ を求めるか？
y(t )データモデル：=yS (t ) + ε において：Q

信号ベクトル yS (t ) = ∑ sq (t )lq の時系列データに注目する．
q =1

種々の時間点でのyS (t1 ), yS (t2 ), , yS (tK ) はK が十分大きければ信号部分空間を埋め尽くす．つまり，以下が成立．
信号部分空間＝ span {l1, , lQ } = span { yS (t1 ), , yS (tK )}=y(t ) yS (t ) + εとは言ってもこれも未知量である．

これのみが知っている量
本文書の著作権は株式会社シグナルアナリシスに帰属します．

信号部分空間： span { y (t ), , y (tS1SK)} を推定する．

信号行列：[ yS (t1 ), , yS (tK )] の特異値展開

γ 1
γQ yS (t1 ), , yS (tK ) =  u1, , uQ , uQ +1,  uM  0T  v1     vT   TQ 0 vQ +1    0   vT  M 0span { yS (t1 ), , yS (tK )} = span {u1, , uQ } であり，u1, , uQ が信号部分空間の正規直交基底である．
しかし， yS (t1 ), , yS (tK ) が未知量なため， u1, , uQも未知量である．

本文書の著作権は株式会社シグナルアナリシスに帰属します．

データモデル=y(t ) yS (t ) + ε
γ 2=R RS + Rε
 12γ
Q
RS =  u1, , uQ , uQ +1, , uM  0

サンプル共分散行列  uT1     uT  Q   uQT +1 0   0   uT  M 
02=ρI であるので，εデータサンプル数無限大の極限で R

γ 2 + ρ 2
 1



22γ+ρ
Q
R =  u1, , uQ , uQ +1, , uM  0
固有ベクトルは同じになる！!

  uT  1     uT  Q   uT   Q +1   ρ 2   uTM 0

ρ2

本文書の著作権は株式会社シグナルアナリシスに帰属します．

データサンプル数無限大の極限で，データ共分散行列 R の Q 個の最大固有値に対応した固有ベクトル u1, , uQ は信号部分空間の基底ベクトルに等しい．
データサンプル数有限の場合，データ共分散行列 R の Q 個の最大固有値に対応した固有ベクトル u1, , uQ は信号部分空間の基底ベクトルの最尤推定解でとなる．
計測データ行列： [ y(t1 ), , y(tK )] を特異値展開し，Q 個の最大特異値に対応したQ 個の空間特異値ベクトルが信号部分空間の正規直交基底の最尤推定解である．

ただし，Q は未知量であり，通常，特異値（固有値）スペクトルから推定する．
本文書の著作権は株式会社シグナルアナリシスに帰属します．

Signal space projection
センサーノイズの除去
計測データ=y(t ) yS (t ) + ε

データ行列 [ y(t1 ), , y(tK )] の特異値展開によりsignal subspaceの基底ベクトルを推定し，プロジェクターP を構成する．
このとき，信号数Q は未知量であり，特異値（固有値）スペクトルのthresholdingにより推定する．
計測データにプロジェクターを乗ずることにより，すなわち，Py(t )= PyS (t ) + Pε= yS (t ) + Pε からノイズ除去を行う．
本文書の著作権は株式会社シグナルアナリシスに帰属します．

Signal space projection
センサーノイズの除去より詳しく見てみると，

γ 1
データ行列のSVD：



γQ

 y(t1 ), , y(tK ) = u1, , uQ , uQ +1,  uM 



B

0MT  v1     vQT  T  vQ +1   γ M   vT  M 0γQ +1= ∑ γ j uj vTj
j =1

計測データにプロジェクターを乗ずると，

PBQT
T=uu
uγuv
γuv
[u,,][,,]∑ j j j∑ j j j11
Q
Q=j 1 =j 1
TM

結局，データ行列を特異値展開して小さな特異値をゼロとする方法--SVDフィルターと呼ばれている--と同じになる
本文書の著作権は株式会社シグナルアナリシスに帰属します．

Signal space projection
妨害信号の除去
データモデル
計測データ

y(t ) = yS (t ) + yI (t ) + ε

妨害磁場データy=(t ) yI (t ) + ε
C

妨害磁場データ空間（interference subspace)を推定するために妨害磁場のみのデータ yC (t ) を必要とする．empty-room noise dataが用いられる．

データ行列 [ yC (t1 ), , yC (tK )] の特異値展開によりinterference subspaceの基底ベクトルを推定し，プロジェクターP を構成する．

本文書の著作権は株式会社シグナルアナリシスに帰属します．

プロジェクターをデータベクトルに作用させると PyI (t ) = yI (t ) であるので( I − P ) y(t ) = ( I − P ) yS (t ) + ( I − P ) yI (t ) = yS (t ) − PyS (t )信号成分の中で妨害磁場空間に含まれる成分最後の項 PyS (t ) が，どれくらい小さくなるかは，妨害磁場空間と信号磁場空間の空間的な分離の程度による．あまり分離していなければ，この操作により信号も小さくなってしまう．
空間的な分離が大きければ，信号の損失は小さくてすむ．

本文書の著作権は株式会社シグナルアナリシスに帰属します．信号部分空間を時間領域で定義する測定が行われた時間点：
q 番目のソースのタイムコースを要素とする行ベクトル：

を時間領域の信号部分空間と定義

信号行列：

各時間点の信号ベクトルを列として並べて作る行列

空間と時間の信号部分空間に対称な関係がある
Each column ofColumn span of B = E
SSEach row of

Row span of B = K
S
S
本文書の著作権は株式会社シグナルアナリシスに帰属します．

時間・空間信号部分空間• 信号部分空間は（従来）信号源リードフィールドベクトルのスパンとして（空間領域で）定義される．
• 時間領域の信号部分空間を信号源のタイムコースを要素とする行ベクトルのスパンで定義する．
• 空間領域の信号部分空間と時間領域の信号部分空間に対照的な関係を見出すことができる．

表記法
行列Xの列空間(column space)：csp(X)
行列Xの列空間(row space)：rsp(X)csp( BS ) = ESrsp( BS ) = K S

本文書の著作権は株式会社シグナルアナリシスに帰属します．

Signal space projection (SSP)(まとめ)Data model: B = B

S

空間領域SSP+ BI + Bε

• 別に取られた妨害磁場データを用いて，妨害磁場部分空間の基底ベクトルを求め，妨害磁場部分へのプロジェクター P を導出．
I

• 妨害磁場除去はデータを妨害磁場空間の直交方向へ投影: ( I − P )BIすることにより行う．
時間領域 SSP• 時間領域の妨害信号部分空間へのプロジェクター除去は B( I − Π ) として行う．

Π I が求まれば，妨害信号I•

すでに提案されているいろいろな方法がタイムドメインSSPと解釈でき，プロジェクターの求め方がそれぞれの方法の違いとなっている．

いろいろな Π I の求め方を見ていくと本文書の著作権は株式会社シグナルアナリシスに帰属します．

Adaptive noise cancelling（計測センサー・リファレンスセンサー１個づつの場合）
２つのセンサーからの出力が以下のモデルに従うと考える(t )
計測センサー： y=s(t ) + w(t )参照センサー： x (t ) = aw (t )：関心対称の信号：妨害信号

未知のゲインa が既知であれば話は簡単，しかし，未知の場合どうすればよいか？

y(t ) α x (t ) + v(t )y(t)をx(t)で回帰する．つまり， =として回帰残差 v(t) が妨害信号除去後の信号である．

ここで，は ∑ (y(t ) − α x (t ))2 を最小にするとして求める．

t

本文書の著作権は株式会社シグナルアナリシスに帰属します．

2∑ (y(t ) − α x (t )) を最小にするは

α = ∑ y(t )x (t ) / ∑ x 2 (t )tttとして求まり,妨害信号除去後の信号は以下で求まる：∑ y(t )x (t )(t ) y(t ) − tv=2∑ x (t )x (t )t

複数個のセンサーの場合
計測センサー： y(t )=y(t ) Ax(t ) + v(t ) として，以下で求まる．

参照センサー： x(t )v=(t ) y(t ) − yx xx−1 x(t )ここで yx=y(t ) xT (t ) xx ∑ x(t ) xT (t )∑=tt

本文書の著作権は株式会社シグナルアナリシスに帰属します．

補足説明：の導出

y(t ) Ax(t ) + v(t )y(t ) を x(t ) で回帰する: =ˆA2

arg min ∑ [ y(t ) − Ax(t )]AtからよりAの最適推定解を求める．

2∂∂
TT TT T∑ [ y(t ) − Ax(t )] =∑ [ − y Ax − x A y + x A Ax ]
∂A t
∂A t
T
T=−2∑ yxT + 2A∑ xxT =0 より∑ [ −2 yx + 2 Axx ] =tt

ˆ ∑==AyxT (∑ xxT )−1 Σyx Σ xx−1t

したがって，ttを得る．
である．

本文書の著作権は株式会社シグナルアナリシスに帰属します．

Adaptive noise cancelling (ANC)（複数センサーの場合）
センサーデータ：リファレンスデータ

y(t )信号をリファレンスセンサー データで回帰する: =Ax(t ) + v(t )妨害信号除去結果，つまり回帰残差は時系列データ行列を導入：

X = [ x(t1 ), , x(tK )]Y = [ y(t1 ), , y(tK )]V = [v(t1 ), , v(tK )]Σyx

1 K11 K
1TTTT()(),()()=ytxt=YXΣ=xtxtXX∑ k∑kxxk
K k 1=K
Kk1 kK
本文書の著作権は株式会社シグナルアナリシスに帰属します．

Adaptive noise cancelling (ANC)
は行列で表現すると，回帰残差:と表される

Xの行空間へのプロジェクター• ANCは計測データをリファレンスセンサーデータ行列の行空間の直交方向に投影する方法である．
• したがって，この行空間が妨害信号部分空間を近似していると考えれば，ANCは時間領域 SSPと解釈できる．
• リファレンスセンサーデータ行列 の行空間は，low-rank signalの仮定のもとで，高SN比であるなら妨害信号部分空間を近似している．
本文書の著作権は株式会社シグナルアナリシスに帰属します．

妨害磁場ソース2個による実験

256-channel CTFセンサーアレイと架空のリファレンスセンサー（６チャンネル）信号磁場センサータイムコース(上段）と妨害磁場の重畳したタイムコース(下段）．妨害磁場は２個のランダムな活動を仮定．

リファレンスセンサーのタイムコース(上段）とANCによる妨害信号除去結果．

リファレンスセンサーのタイムコースに低周波の周期ノイズを混入．

リファレンスセンサーに混入しているノイズは結果に影響しない．
本文書の著作権は株式会社シグナルアナリシスに帰属します．

妨害磁場ソース６個による実験

妨害磁場の重畳したタイムコース．妨害磁場は６個のランダムな活動を仮定して発生．

リファレンスセンサーのタイムコース(上段）．リファレン
スセンサーのタイムコースに低周波の周期ノイズを混入．
ANCによる妨害信号除去結果（下段）リファレンスセンサーのタイムコース(上段）．リファレンス
センサーのタイムコースに低周波の周期ノイズが混入していない場合．ANCによる妨害信号除去結果（下段）Low-rank信号の仮定が成り立たない場合にはリファレンスセンサーに混入しているノイズがANCの結果に影響する．
本文書の著作権は株式会社シグナルアナリシスに帰属します．

シミュレーション結果の解釈
推定された
妨害信号部分空間

妨害信号部分空間
信号部分空間M次元ベクトル空間

信号部分空間真の妨害信号部分空間

M次元ベクトル空間

リファレンスセンサー出力に，計測センサー出力には含まれない活動が存在す
ると妨害信号部分空間は真のものより広く推定される．しかし「推定された妨害信号部分空間」の補空間は十分広いため，「真の妨害信号部分空間」の補空間とさほど違わずに信号部分空間を含み，妨害信号除去結果はそれほど影響を受けない．

本文書の著作権は株式会社シグナルアナリシスに帰属します．

シミュレーション結果の解釈:low-rank signal の仮定が当てはまらない場合
妨害信号部分空間
信号部分空間M次元ベクトル空間

信号部分空間真の妨害信号部分空間

推定された
妨害信号部分空間
M次元ベクトル空間

low-rank信号の仮定が当てはまらないと妨害信号部分空間の推定は間違ったものになり，「推定された妨害信号部分空間」の補空間は真の妨害信号部分空間を含むため，ANCを用いても妨害信号は取りきれない．
本文書の著作権は株式会社シグナルアナリシスに帰属します．

Common temporal mode subspace projection (CtSP)
Data modelSensor data：

B = BS + BI + Bε
 = B + B + BReference sensor data： BIwε

これらデータから（時間領域）妨害信号部分空間 K I = rsp( BI )
する．

を推定

センサーノイズ項を無視して，

B BS + BI=B B I + B w
=rsp( B) ⊂ rsp( BS ) + rsp( BI ) =KS + KIrsp( B ) ⊂ rsp( B I ) + rsp( B w ) =+KI Kw
本文書の著作権は株式会社シグナルアナリシスに帰属します．

Common temporal mode subspace projection (CtSP)
若干長い導出の後：KI⊃ rsp( B)  rsp( B ) を得る．

rsp( B)  rsp( B ) の基底ベクトルを求めてプロジェクター Π I を計算し，

ˆ=BB( I − Π I ) として時間領域SSPを行い，妨害信号を除去する．
SK I = rsp( B)  rsp( B ) が成立しないので，理論的には妨害信号を完璧に除去できるわけではない．

本文書の著作権は株式会社シグナルアナリシスに帰属します．

Dual signal subspace projection (DSSP)
Data model

B = BS + BI + Bε

リファレンスセンサーデータを必要としない．なぜなら，計測データ B から，信号を含まない，リファレンスセンサーデータを作り出す．
どうするか？
信号部分空間へのプロジェクター P が既知なら，P BS

S

S= BS であるので，( I − PS )B = ( I − PS )BI + ( I − PS )Bεから，synthetic reference dataを作り出せる．
信号部分空間へのプロジェクターは未知なため，次の疑似信号部分空間プロジェ
クターを代わりに用いる．

本文書の著作権は株式会社シグナルアナリシスに帰属します．

疑似信号部分空間プロジェクター

rj

ボクセルリードフィールド行列
LV =  L(r1 ),, L(rN ) SVD of LVが他の特異値と比べ有意に大きいとすれば，は近似的に信号部分空間を構成する．
T =  u ,, u   u ,, u  により，近似的に，
プロジェクター： Pξ 1ξ 1

 = B が成立．
PBSS
本文書の著作権は株式会社シグナルアナリシスに帰属します．

Dual signal subspace projection (DSSP)
Data model

B = BS + BI + Bε

データから（時間領域）妨害信号部分空間 K I = rsp( BI ) を推定する．

 と I − P を乗じて以下を得る．
データ行列に P
 = PB + PB + PBPBSIε( I − P )B = ( I − P )BI + ( I − P )Bε ) ⊂ rsp( PB ) + rsp( PB ) + + rsp( PB =rsp( PBSIε ) KS + KI + Kεrsp(( I − P )B ) ⊂ rsp(( I − P )B I ) + rsp(( I − P )B ε ) =K I + K ε' )  rsp(( I − P )B)
若干長い導出の後： K I ⊃ rsp( PBを示すことができる．

本文書の著作権は株式会社シグナルアナリシスに帰属します．

Dual signal subspace projection (DSSP) )  rsp(( I − P )B) の基底ベクトルを求め，これらから計算したrsp( PB
プロジェクターを，妨害信号部分空間へのプロジェクター Πˆ
時間領域SSP =BSKIIとして，B( I − Π I ) を実行し，妨害信号を除去する．

 )  rsp(( I − P )B)が成立しないので，理論的にはrsp( PB妨害信号を完璧に除去できるわけではない．

B = BS + BI + Bε において BI が大きく支配的な場合には，
 )  rsp(( I − P )B) が成り立っている（と思われる）．
K I ≈ rsp( PB

本文書の著作権は株式会社シグナルアナリシスに帰属します．

部分空間FとGの共通部分 F  G の基底ベクトルの求め方
以下に詳細な記載あり

資料：「線形数学の基本事項」の12節「部分空間の考え方を基にした環境ノイズ除去アルゴリズム」のAppendix B
FとGの大きさ１に規格化された要素ベクトルを f1 と g1 とする．
subspace angle:

cos(θ1 ) = max f1T g1
f1∈F ,g1∈G=cos(θ2 )T
to f2T f1 0,=
and g 2T g1 0=max f2 g 2 subject

f2∈F ,g 2∈Gとして繰り返すことにより，principal angle:
principal vector:

θ1 ≤ θ2 ≤  ≤ θνf1, f2 ,  , fν g1, g 2 ,  , gνが求まる．principal angle は2つの部分空間の「離れ具合」あるいは「一致さ」を表す． ν は2つの部分空間の小さなほうの次元を表す．
本文書の著作権は株式会社シグナルアナリシスに帰属します．

直感的な例
FとGが角度 φ で交わる直線

直線と平面

θ1 = φθ1 = φ

２つの平面3次元空間と平面=
θ1 0,=
θ2 φ
2つの部分空間の共通集合を見つけるにはこれら部分空間のprincipal angleを計算する．

θ1 = =θr =0 ≤ θr +1 ≤  ≤ θνを見出せば共通集合の次元はrであり，その正規直交基底ベクトルは=
θ1 0,=
θ2 0

f1, f2 ,  , fr あるいは g1, g 2 ,  , g r で与えられる．
本文書の著作権は株式会社シグナルアナリシスに帰属します．

どうやって，共通集合の正規直交基底を求めるか部分集合FとGの正規直交規定を行列：U = [ u ,  , u ]µ1u1,  , uµ ，v1,  , vν とする．

V = [v1,  , vν ] とすれば，FとGの任意の要素ベクトルを f と g は以下のように表される．
µc j uj Uc, c [c1, , c µ ]T
f = ∑==j =1

νg= ∑==d j v j Vd, d [d1, , dν ]Tj =1cos(θ1 ) = max f T gf ∈F ,g ∈GTTcos(θ1 ) = max cT [U TV ]dc,dmax c [U V ]d は特異値の性質から，行列U V の最大特異値に等しい．
Tc,dprincipal vector f1, f2 ,  , fr は，r 個の最大特異値に対応した特異値ベクトルc j ( j = 1, , r ) を用いて f j = Uc j ( j = 1, , r ) で与えられる．
本文書の著作権は株式会社シグナルアナリシスに帰属します．

データ行列X とY に対しrsp(X)とrsp(Y)の共通集合の基底ベクトルを求める方法

X とY を特異値分解し，行方向の特異値ベクトルから，行空間の正規直交基底を求める．

rsp( X ) = rsp([ x1T , , xTµ ]T )rsp(Y ) = rsp([ y1T , , yνT ]T )U[=x1T , , xTµ ], V [ y1T , , yνT ] (µ > ν ) としてU とV の列空間の主角度は，行列：U TV の特異値で与えられる．特異値が以下の特性：

γ 1 =  = γ r ≈ 1 > γ r +1 ≥  ≥ γ νを持てば， rsp (X)とrsp(Y)の共通集合の次元を r と決めることができる．
共通集合の基底

ψ1, , ψr は U TV の特異値ベクトルを uj , v j として

ψk = Uuj あるいは=
ψk Vv=( j 1, , r ) で与えられる．
j
本文書の著作権は株式会社シグナルアナリシスに帰属します．

Data modelSensor noise suppression (SNS)Sensor data： =B

BS + BεSNS法の仮定：
あるセンサーのタイムコースは他のセンサーのタイムコースの線形和で表される．

データ行列B の j 番目の行を β j と書く．
行列： Θj = [ β1 , , β j −1, β j +1, , βM ]T

重みをTTTと定義し，wi と書けば，T,[,,,,,]=βTj=wi βiT Θ=www
ww
w
∑j
1j −1j +1Mi≠jである．
本文書の著作権は株式会社シグナルアナリシスに帰属します．

重みベクトル w は最小二乗法から

w = (ΘTj Θj )−1 ΘTj βTjと求まるので，SNS法による j 番目のセンサーのノイズ除去結果は

βˆj = β j Θj (ΘTj Θj )−1 ΘTjとして求まる．つまり，オリジナルタイムコース β を rsp(Θj ) に投影するこjとによりノイズ除去を行う．
であるのでrsp(Θj ) ≈ rsp( B) ≈ rsp( BS ) =KS
SN比の高い場合
SNSアルゴリズムは信号部分空間を rsp(Θj ) で近似した時間領域SSPであると解釈できる．
本文書の著作権は株式会社シグナルアナリシスに帰属します．