圧縮センシングの遠隔制御系への応用
Compressed Sensing Application to Remote Control Systems
永原正章 †松田崇弘 ‡林和則 †† 京都大学大学院情報学研究科‡ 大阪大学大学院工学研究科

Masaaki Nagahara†
Takahiro Matsuda‡
Kazunori Hayashi††Graduate School of Informatics, Kyoto University‡Graduate School of Engineering, Osaka University
アブストラクトIn remote controlled systems, eﬃcient
representation of control signals is one of the crucial is-れる場合が多い．また，通信路で雑音等の原因によりビッsues because of the bandwidth-limited communication
channel, such as a wireless communication link, between

信号を送る必要があり，制御入力信号に対する伝送帯域

the controller and the controlled object. In this article, we propose a compressed sensing-based sparse rep-帯域で効果的に制御信号を送るための手法が重要である．

resentation of control signal, which enables us to reducethe diﬀerential entropy of the control signal to be transmitted from the controller to the plant, while keeping

ト誤りが生じる場合には，誤り訂正符号等で冗長化したの制約はより強いものとなる．従って，与えられた伝送
本稿では，圧縮センシングの手法を用いてスパースな（すなわち多くの要素が 0 であるような）制御入力信号を得ることにより，これを表現するのに必要な情報量を低減する手法を提案する．本手法は制御入力信号の情報圧縮

the controlled performance. The sparse representation
can be eﬀectively obtained by solving l1 -l2 optimization

法の一つととらえることができるが，通常の情報源符号

problem with the iterative shrinkage method. We show
design examples to illustrate the eﬀectiveness of the pro-化しないように情報源の確率分布そのものを変更するこ化における情報圧縮とは本質的に異なり，制御性能が劣とで，その微分エントロピーを減少させるものである．

posed method.
1

具体的には，制御対象の出力軌道に対する追従誤差を2

 ノルムで測り，かつ送信信号（設計パラメータ）のスパースネスに対して 1 ノルムのペナルティ項を加えるこ

はじめに圧縮センシングまたは圧縮サンプリングは，元のデーとにより，制御性能の向上とスパースネスとを同時に達タがスパースであるという仮定を置くことによって，非成する制御信号を設計する．この設計は，1 ではなく 2

常に少ない観測データから元のデータをかなり正確に復のペナルティを考えた場合，S. Sun ら [4] によって提案さ元できる手法として，近年，とくに信号処理や情報理論れた制御理論的スプラインの問題となる．一方，我々のの分野で盛んに研究されている [2], [5]．この手法の特徴

問題は，1 /2 最適化として定式化されるが，上記で述べは，スパースネスの指標であるベクトルの非ゼロ要素の数た繰り返し縮小法により効率的に最適解を得ることがで0（しばしば “ ノルム” と呼ばれる）の最小化を，その凸きる．この設計によって，制御性能を大幅に劣化させる緩和である 1 ノルムの最小化に置き換えて逆問題を解くことなく，制御信号の微分エントロピーを減少させるこ

1点にある．これは，最小二乗法の正則化項に  ノルムをとが可能であることを設計例により示す．

用いる LASSO (Least Absolute Shrinkage and SelectionOperator) と呼ばれるものとも関係がある [5]．この最適解は，内点法を用いて計算することもできるが，近年，繰り返し縮小法 (iterative shrinkage) と呼ばれる手法が提案され，より高速に最適解を見つけることが可能となった記法
ベクトル v = [v1 , . . . , vn ] ∈ Rn に対して，その 1 ノルムおよび 2 ノルムを[1], [6]．
一方遠隔制御では，制御入力信号は無線伝送路などの伝送帯域の制約が強い通信路を介してシステムに加えらv1 :=n
i=1|vi |,v2 := n
i=11/2vi2で定義する．また行列 Φ に対して，そのノルムを

D

Φv2
Φ := max= σmax (Φ)
v2 =0 v2で表す．ここで σmax (Φ) は行列 Φ の最大特異値である．
実数 x ∈ R に対して

sgn(x) :=⎧⎨1,⎩−1,(λW −1 + G0 )−1

gi (t)uPy

図 1: Remote control system optimized with J2 (u) in (2).

The vector θ is transmitted through a communicationchannel.

if x ≥ 0,
if x < 0,関数 gi (t) は(x)+ := max{x, 0}
gi (t) :=と定義する．区間 [0, T ] (T > 0) 上の 実数値 Lebesgue

2 乗可積分関数の全体を L2 [0, T ] と書く．f, g ∈ L2 [0, T ]⎧⎨CeA(ti −t) B,

if ti > t,⎩0,

if ti ≤ t(4)で定義され，係数 θ := [θ1 , . . . , θN ] はに対して，内積をf, g := T

θ = (λ2 W −1 + G0 )−1 αf (t)g(t)dt0

問題設定
次の線形システム P を考える:
ẋ(t) = Ax(t) + Bu(t),P :y(t) = Cx(t),(5)で与えられる．ここで，α := [α1 , . . . , αN ] であり，また⎡⎤gt1 , gt1  . . . gtN , gt1 ⎢
⎥
..
..
..
⎥,G0 := ⎢.
.
.
⎣⎦gt1 , gtN  . . . gtN , gtN で表す．

2

θt ∈ [0, ∞).

W := diag {w1 , . . . , wN }(1)である．遠隔制御系の枠組みで考えると，この制御では，
ここで A ∈ Rn×n , B ∈ Rn×1 and C ∈ R1×n とし， まずデータ D が与えられたとき，(5) によって係数ベクx(0) = 0 とする．また，(A, B) は可制御，(C, A) は可 トル θ が計算され，この θ が通信路に送信される．受信観測とする．このシステムの出力信号 y に対して，次のデータ量は少ない方が望ましい．したがって，圧縮セン

D := {(t1 , α1 ), (t2 , α1 ), . . . (tN , αN )}.

シングの考え方 [2], [5] に従い，送信信号 θ のスパース表
現を考える．まず，制御入力 u は，L2 [0, T ] の部分空間
M
2
VM := u ∈ L [0, T ] : u =θ i ψi , θ i ∈ R .
(6)

ここで，各 ti は時刻を表すデータで0 < t1 < t2 < · · · < tN =: Tとし，その各時刻において，システムの出力 y(ti ) が αiに近くなるように，すなわち y(ti ) ≈ αi となるように制御入力 u を設計したい．
この設計問題に対して，次の 2 ノルム評価関数を考える:

J2 (u) := λ2側では，(3) により制御入力信号 u が計算される．この状況では，通信路の容量に制約があるのが通常で，θ のデータ点を与える:

 Tu(t)2 dt +0Ni=1に属する信号と仮定する（M ≥ N とする）．ここで，{ψi }は L2 [0, T ] の線形独立なベクトルとする．このとき，評価関数 J2 (u) の第二項は次のように表現できる:N

wi (y(ti ) − αi )2 .

(2)i=1i=1= θ G W Gθ − 2θ G W α + α W α

ここで，wi > 0 は各データ点へのフィッティングに対す

2= Φθ − v2 .

る重みである．この評価関数 J2 (u) を最小化する最適制
御入力は，次式で表されることが知られており，制御理論的スプラインと呼ばれる [3], [4]:u(t) =N

i=1

θi gi (t).

(3)

wi (y(ti ) − αi )2

ただし，

⎡

g1 , ψ1 ⎢..
⎢G := ⎣.
g1 , ψM ...
..
.
...

⎤
gN , ψ1 ⎥..
⎥.
⎦gN , ψM 

とし，Φ := W 1/2 G,v := W 1/2 α とする．これより，係θ数ベクトル θ に 1 のペナルティを与えた次の評価関数を導入する:
2J1 (θ) := λθ1 + Φθ − v2 .

(7)

−λ/c

上で述べたように，第一項の 1 ノルムはベクトル θ のスη
0λ/c

パースネスを表す指標である．この評価関数を最小化する最適解を θ∗ とすると，最適制御入力はu∗ (t) =Nθi∗ ψi (t),t ∈ [0, T ]i=1で与えられる．

3

図 2: Nonlinear function θ = sgn(η)(|η| − λ/c)+

1 -2 最適化によるスパース表現(7) の評価関数は θ に関して凸であり，最適解は必ず存在する．しかし，その最適解は (5) のように解析的にはlimk→∞ θ[k] は (7) を最小化する最適解に一致する [1]．

表すことができず，最適解を求めるためには繰り返し計

さらに，最適解への収束の早さは最悪評価で O(1/k 2 ) で算が必要である．この繰り返し計算には，内点法などがあることも知られている [1]．上記のアルゴリズムは，計従来用いられてきたが，近年，繰り返し縮小法 (iterative算機に容易に実装が可能であり，収束も早いことから，リshrinkage)[1], [6] と呼ばれる手法が提案され，これによりアルタイムでスパースなベクトル θ を得ることが可能で高速な数値計算が可能となった．

ある．

繰り返し縮小法のアルゴリズムは以下で与えられる [1],[6]: 初期値 θ[0] ∈ RN を与え，β[1] = 1, η[1] = θ[0] とする．定数 c を4設計例

c > Φ2 = σmax (Φ)2となるように選ぶ．次を繰り返す1 :
1 
Φ (v − Φη[k]) + η[k]
θ[k] = Sλ/cc
1 + 1 + 4β[k]2
β[k + 1] =2
β[k] − 1η[k + 1] = θ[k] +(θ[k] − θ[k − 1])β[k + 1]

k = 1, 2, . . .

⎢Sλ/c (η) := ⎢⎣⎡01000⎢A = ⎣0⎤0⎥
1⎦ ,0⎡ ⎤0⎢ ⎥B = ⎣0⎦ ,1
C= 10


0 .とする．すなわち，P (s) = 1/s3 である．データ D のサンプル時刻を ti = 0.1 + (i − 1) × 0.5, i = 1, 2, . . . , 11}関数 Sλ/c は η = [η1 , . . . , ηN ] に対して⎡

線形システム P の実現を

⎤sgn(η1 )(|η1 | − λ/c)+⎥..
⎥.
⎦sgn(ηN )(|ηN | − λ/c)+とし，データ αi を正弦波 αi = sin ti とする（図 3 を参照）．また，L2 [0, T ] の線形独立なベクトル ψi (t) は，上記の A, B, C および {t1 , . . . , t11 } を使って (4) によって定義される gi (t) と同じものとする2 ．また，(2) に対するで定義される非線形関数である．関数 Sλ/c の各要素に対λ を 0.1，(7) に対する λ を 0.02 とする．これらの設定のもとで，(3) で計算される， (2) の 2 評価関数を最適する非線形関数 sgn(η)(|η| − λ/c) の入出力関係を Fig. 2化する θ を θ2∗ とおき，また (7) の 1 -2 評価関数を最適2に示す．条件 c > Φ が成り立つような c > 0 を選べば，上のアルゴリズムは任意の Φ と v に対して収束し，1 繰り返し縮小法にも様々な方法が提案されているが [6]，ここで用いた方法は FISTA (Fast Iterative Shrinkage-Thresholding Algorithm)と呼ばれる方法である [1]．

2 (A, B) が可制御かつ (C, A) が可観測ならば，(4) で定義される{gi } は線形独立となることが知られている [3]．

Original curve and sampled data
10.80.6
y0.4
10.2
00.5
-0.2
-0.4
0-0.6
-0.8-0.5-1

0

12

3

456-1

図 3: Original curve sin t and sampled data D.

-1.50

1

23456

図 4: Original curve sin t: (blue line), Output y(t) of化する θ を θ1∗ とおくと，それぞれ以下で与えられる．

⎤0.9919⎥⎢⎢ 4.7303 ⎥⎥⎢⎢ 5.4459 ⎥⎥⎢⎥⎢
⎢ 3.9878 ⎥⎥⎢⎢ 1.8978 ⎥⎥⎢⎥⎢θ2∗ = ⎢ 0.1898 ⎥ ,⎥⎢⎢−0.8228⎥⎥⎢⎥⎢⎢−1.2059⎥⎥⎢⎢−1.0451⎥⎥⎢⎥⎢−0.3015⎦⎣
1.0642⎡⎡00

1 -2 optimal input with λ = 0.02 (red line), 2 optimal
input with λ = 0.1 (blue dots).

⎤⎥
⎢
⎥
⎢
⎥
⎢
⎢ 2.0665 ⎥⎥⎢⎥⎢⎢ 3.8879 ⎥⎥⎢⎢ 3.6590 ⎥⎥⎢⎥⎢θ1∗ = ⎢ 1.1288 ⎥⎥⎢
⎥
⎢
0⎥
⎢
⎥
⎢
⎢−1.7358⎥⎥⎢⎢−2.6657⎥⎥⎢⎥⎢−0.3231⎦⎣1.7432

u543210-1

1 -2 最適化によりスパースな解 θ1 が得られることが-20123456わかる．これらのパラメータに対するシステムの出力をFig. 4 に，入力 u(t) を Fig. 5 に示す．また復元誤差
e(t) = | sin t − y(t)| を Fig. 6 に示す．1 -2 による最適化
2によるスパース表現にもかかわらず， 最適化とほぼ同じ性能を示していることがわかる．

さらに，(7) の 1 -2 評価関数の λ を 0.05 にした場合

図 5: 1 -2 optimal input u(t) with λ = 0.02 (red line),2 optimal input u(t) with λ = 0.1 (blue dots).

error
0.70.60.5u50.4
40.3
30.2
20.1

0

101234560-1

図 6: Absolute error with θ1∗ (solid) and θ2∗ (dots).の最適パラメータ θ1∗ は以下のようになる．
⎤
⎡
0⎥
⎢
⎥
0⎢
⎥
⎢
⎢ 0.3672 ⎥⎥⎢⎥⎢⎢ 3.2113 ⎥⎥⎢⎢ 3.2282 ⎥⎥⎢⎥⎢θ1∗ = ⎢
⎥.
0⎥
⎢
⎥
⎢
0⎥
⎢
⎥
⎢
⎢−0.6546⎥⎥⎢⎢−2.2377⎥
⎥
⎢
⎥
⎢
0
⎦
⎣

-2

0

1

2

3456

図 7: 1 -2 optimal input u(t) with λ = 0.05 (red line),2 optimal input u(t) with λ = 0.1 (blue dots).

1.0940error

定数 λ の値を大きくしたことで，よりスパースな表現が得られていることがわかる．この θ1∗ を用いた入力 u(t)，および復元誤差をそれぞれ Fig. 7，Fig. 6 に示す． 入力u(t) は λ = 0.02 の場合に比べて小さくなっていることがわかる．これは，1 ノルムのペナルティをより大きくとったことに対応している．さらに，このペナルティを大0.7

0.6

0.5

0.40.3きくすることにより，復元誤差が上昇していることがわ
かる．この例題から，入力信号のパラメータ θ のスパー
スネスと復元誤差（性能）との間にはトレードオフの関0.20.1係があることがわかる．
次に信号 θ のスパース性と情報圧縮との関係を調べる

00123456

ために，(1) のシステムをランダムに 1000 サンプル発生させる．ここで，システムの次数は n = 2 に固定する．なお，システムに対するランダムサンプルの発生に2は SCILAB の ssrand 関数を用いた． 最適なベクトル θ2∗ と 1 -2 最適なベクトル θ1∗ の要素の分布を調べる．
得られた θ の確率分布を図 9 に示す．なお，1 -2 最適化によって得られた分布は原点 θ = 0 で急峻に立ち上が図 8: Absolute error with θ1∗ (solid) and θ2∗ (dots).

参考文献

Probability Distribution
0[1] A. Beck and M. Teboulle, A fast iterative shrinkage-−0.5

thresholding algorithm for linear inverse problems,
SIAM J. Imaging Sciences, vol. 2, No. 1, 183/202−1P(θ)−1.5(2009)log

10−2−2.5[2] E. J. Cnades and M. B. Wakin, An introduction
to compressive sampling, IEEE Signal Processing
Magazine, vo. 25, 21/30 (2008)−3−3.5−4−4.5

−0.8

−0.6

−0.4

−0.20θ

0.20.40.60.8

図 9: Probability distribution log10 P (θ): 2 optimal(dots) and 1 -2 optimal (solid). A logarithmic scale
is used.

る特性（デルタ関数のような特性）を持つため，比較を[3] M. Egerstedt and C.F. Martin, Control Theoretic
Splines, Princeton University Press (2010)[4] S. Sun, M. B. Egerstedt, and C. F. Martin, Control
theoretic smoothing splines, IEEE Trans. on Automatics Control, Vol. 45, No. 12, 2271/2279 (2000)[5] 田中，圧縮センシングの数理，IEICE Fundamental
Review, Vol. 4, No. 1, 39/47 (2010)

容易にするために log10 P (θ) を図示している．これより，[6] M. Zibulevsky and M. Elad, L1-L2 optimization in

集中しており，微分エントロピーが小さくなること予想

signal and image processing, IEEE Signal Processing Magazine, Vol. 27, 76/88 (2010)

2 最適化によって得られた分布よりも原点付近に分布がされる．実際，この分布にもとづき微分エントロピーを計算すると，

H1 (θ) = 2.4885,H2 (θ) = 7.3817となる．ただし，H1 は 1 -2 最適化によって得られた θの微分エントロピー，H2 は 2 最適化によって得られた

θ の微分エントロピーである．これより，1 -2 最適化による微分エントロピーが 2 最適化の微分エントロピーの約 1/3 になっていることがわかる．したがって，1 -2 最適化によるスパース表現により，送信信号のより大きな情報圧縮が可能となる．

5

おわりに
本稿では，遠隔制御系の枠組みで，送信信号（最適制御入力を定義するパラメータ θ）のデータ量を抑えるために，1 -2 最適化の枠組みで問題を定式化し，繰り返し縮小法により効率的に最適解を求めることができることを示した．また設計例を示し，性能をそれほど下げずにある程度のスパース表現が可能なこと，またパラメータのスパースネスと制御性能の間にトレードオフの関係があることを示した．また，スパース表現により，送信信号の微分エントロピーが小さくなることを数値例によって示し，送信信号の大きな情報圧縮が可能であることを示した．