はじめに

こんにちは。(株) 日立製作所の Lumada Data Science Lab. の松本茂紀です。

IoTプラレールを作ってみましたシリーズの第5弾です。

今回は「プラレールとフィットネスバイクを電波でつなごう」です。

その他の記事は以下から参照できますのでよろしければどうぞ。

IoTプラレールを作ってみました（その1：IoT電池の作成前編）

IoTプラレールを作ってみました（その2：IoT電池の作成後編）

IoTプラレールを作ってみました（その3：IoTサイクルメータを作ろう前編）

IoTプラレールを作ってみました（その4：IoTサイクルメータを作ろう後編）

IoTプラレールを作ってみました（その6：プラレールとフィットネスバイクを電波でつなごう後編）

要件定義

今回の記事では、前回の記事で登場したパルスオキシメータの残りの実装と、いよいよプラレールとフィットネスバイクとの間で通信させるお話です。今回もまずは重要な要件定義から始めます。これまでの記事でも既に触れた２つの要件を満たすような機能を検討します。

運動による消費カロリーをおおよそ知りたい
プラレールで長い時間遊びたい

１つ目の要件は、前回の記事で作成した指サック型のパルスオキシメータを使って心拍を測定し、そこから消費カロリーを算出したいというものです。正確に消費カロリーを知るのは、本来は非常に大変は計測が必要ですが、なるべく手間を掛けずに尤もらしい結果が欲しいところです。また、１日の運動時間の目安にしたいので、リアルタイムに何カロリー消費したかの計算ができると嬉しいです。

２つ目の要件は、電池交換が大変手間なので、なるべく長い時間遊べるようにしたいというものです。無線通信によって電池の消費量が著しく低下するのは避けたい所です。そこで通信方式についても検討の必要がありそうです。

パルスオキシメータの原理を知ろう

ここでの目標は、消費カロリーを計算することです。そこに至るまでに必要な、以下のことについて説明していきます。

そもそもパルスオキシメータって何？
心拍数はどうやって計算する？
消費カロリーはどうやって計算する？

そもそもパルスオキシメータとは、動脈血の酸素飽和度(SpO2)や脈拍数を測定する装置で、光を使った血流の変化を調べる方法（Photoplethysmography: PPG）の一つです。パルスオキシメータでは赤色の可視光や赤外線(IR)が利用されています。

赤色光やIRを使ってどうやって酸素飽和度や脈拍数が測れるのでしょうか？その答えは、血液の色に関係があります。動脈の血液は鮮やかな赤色をしていますが、静脈では黒っぽい赤色をしています。赤血球に含まれるヘモグロビンは、ヘム（ヘム鉄錯体）とグロビンから構成されていますが、このヘムが酸素と結びつくことで赤く見えています。

物質が赤く見えるというのは、赤色の可視光だけが物質に吸収されず反射されていることを意味しています。白色に見える光には様々な波長の光が含まれており、光が物質にあたり反射した光の成分によって、その物質の色を認識しています。一般に、目に見える色は、その色の補色にあたる光を物質が吸収していると考えられます。ちなみに、物質が白く見える場合はほとんど光を吸収せず、逆に黒の場合は全ての光を吸収している状態です。私の専門分野の一つである材料科学に関する研究では、物質が吸収する色の違いを、量子力学に基づく分子レベルのシミュレーションによって予測していたため、パルスオキシメータの原理を理解する上で知識が生かされました（笑）

酸素と結合したヘムは、光の波長に対する吸収度合い（吸光度）が600～700 nm付近で著しく低いことが知られいます。まさに、赤色の可視光が吸収されにくく、反射されやすいということです。なお、IR(870~900nm)は酸素の結合有無に対する変化は赤色光に比べ小さいものの、酸素が結合した場合のほうが吸光度が高くなるという、赤色光と逆の傾向が見られます（参考）。

今回用いたパルスオキシメータは、図のように赤色光や赤外光を体に照射し、反射する光の強度を測定しています。光を吸収する体細胞の成分はいくつかありますが、動脈血の層は脈打つたびに血流量が増加し、吸収される光量も増えるため、反射光の強度が減ります。この交流成分を調べることで、心拍数を計算することができます。

なお、ヘモグロビンは４つのヘムを持っており、結合した酸素の数だけ赤色光が吸収されるので、酸素飽和度（ヘモグロビンが結合できる酸素量の何％が結合しているか）を知ることができるというわけです。ただし、酸素飽和度を正確に計算しようとすると、補正など考慮することが多いので、別の機会に検討しようと思います。

パルスオキシメータを使ってみよう

まずはセンサからどんなデータが取得できるかを見てみます。使用したマキシム社のMAX30102センサの情報は、データシートに書いてあります（リンク）。自分で電子工作する場合、モジュールごとのデータシートを確認する必要がありますが、慣れないと読み解くのが大変です。ここでは、ポイントを簡単に説明していきます。

センサからデータを取得するには、レジスタマップへのアクセスが必要になります。MAX30102のレジスタマップは主に４つのカテゴリ（SATUS, CONFIGURATION, DIE TEMPERATURE, PART ID）に分類されています。DIE TEMPERATUREは、LEDから出力される光の波長に影響を与える基板温度に関するもので、SpO2を精度良く算出する際には重要になる項目です。今回は心拍数のみ測るため、必要な項目はSTATUSとCONFIGURATIONだけです。

STATUSとCONFIGURATIONの内容を簡単に説明します。

STATUS

Interruptの情報
センサデータを格納するレジスタに、新たにデータが格納されたときや、データがいっぱいになった時などに、前回記事で接続したGPIO 17番ピンの信号を1にします。このタイミングをトリガーに、データを読み取る処理ができます。今回は一定間隔でデータを読み出すのでこの項目は不要です。
FIFO(First In First Out)の情報
実際のセンサデータが格納されるレジスタ関連の情報です。1回のセンシングで1つのLEDから3バイトの情報が順番に格納されます。最大で32サンプル格納でき、格納されているサンプル数はRead PointerからWrite Pointerの差として得られます。レジスタに格納された順番でデータを読み取る事ができます(先入れ先出し)。

CONFIGURATION

FIFOの設定
格納データのサンプル平均数の設定や、データがいっぱいになったときのアラートの設定ができます。サンプル平均数は４に設定しました。アラート設定は、今回interruptを使用しないので、最大の32サンプルに設定しました。
動作モードの設定
シャットダウン、リセットや、心拍数のみか心拍数とSpO2両方測定するか、などの設定ができます。動作モードを指定する際には、心拍数のみを指定しました。なお、データシートには心拍数のみの場合、アクティブなのは赤色LEDのみと書かれていますが、実際にはIR LEDが動作しています。これに気づくまでは、赤く光らず故障かな？と焦りました。
アナログ-デジタル変換(ADC)の設定
センサの読み取り範囲、サンプリングレート、パルス幅、読み取り解像度を指定できます。読み取り範囲はフルスケールで4096 nA、サンプリングレートは100 Hz、ADC解像度は18ビットに設定しました。なお、FIFOの設定でサンプル平均数を4としているので、レジスタに格納されるデータ数は毎秒25個となります。
LED出力の設定
LEDの電流値を設定します。電流値により、光の強度や、消費電力、温度変化などが変わってくるので、使用条件に合わせて設定します。ここでは4.8 mAに設定しました。

上記をもとに、まずは生のデータを取得するプログラムは以下の通りです。後で、本体プログラムにモジュールとして取り込むことを想定しています。なお、レジスタの書き込み、読み込みの値は、データシートのレジスタマップのビット列表と比較しやすいように２進数で表記しています。このファイルをmax30102.pyという名前で保存します。

import smbus
import time

class MAX30102():

    def __init__(self,channel=1,addr=0x57):
        # デフォルト設定
        #   channel = 1 (Raspberry Piではchannelは1の模様)
        #   アドレスは0x57  (sudo i2cdetect -y 1で確認）

        self.i2c = smbus.SMBus(channel) 
        self.addr = addr 

        #アドレス定義
        self.fifo_w_ptr = 0x04
        self.fifo_r_ptr = 0x06
        self.fifo_data = 0x07
        self.fifo_cfg = 0x08
        self.mode_cfg = 0x09
        self.spo2_cfg = 0x0a
        self.led_1_amp = 0x0c
        self.led_2_amp = 0x0d

        self.configure()

    def configure(self):
        #1.設定をリセット
        self.i2c.write_byte_data(self.addr, self.mode_cfg, 0b01000000 )  #reset

        #2. FIFOのポインタリセット
        self.i2c.write_byte_data(self.addr, self.fifo_w_ptr, 0b00000000 )
        self.i2c.write_byte_data(self.addr, self.fifo_r_ptr, 0b00000000 )

        #3. FIFOの設定：
        #     サンプル平均=4
        #     トリガーは今回使わない
        self.i2c.write_byte_data(self.addr, self.fifo_cfg, 0b01000000 )

        #4. MODEの設定：Heart Rate mode (IR only)
        self.i2c.write_byte_data(self.addr, self.mode_cfg, 0b00000010 )

        #5. SpO2の設定:
        #      ADC range=15.63uA, 4096nA
        #      Sample Rate=100Hz #サンプル平均が４なので 100/4 => 25 samples/sec
        #      Pulse Width=411us
        self.i2c.write_byte_data(self.addr, self.spo2_cfg, 0b00100111 )

        #6. LED Pulse Amplitudeの設定: 4.8 mA
        self.i2c.write_byte_data(self.addr, self.led_1_amp, 0b00100100 )
        self.i2c.write_byte_data(self.addr, self.led_2_amp, 0b00100100 )

    def read_fifo(self):
        # FIFOデータの読み取り
        w_ptr = self.i2c.read_byte_data(self.addr, self.fifo_w_ptr)
        r_ptr = self.i2c.read_byte_data(self.addr, self.fifo_r_ptr)
        num_samples = (w_ptr - r_ptr)%32 #wrap around
        ir_led = []
        while num_samples > 0:
            num_samples-=1
            data = self.i2c.read_i2c_block_data(self.addr, self.fifo_data,3)
            ir_led.append((data[0]<<16 | data[1]<<8 | data[2]) & 0b000000111111111111111111)
        return(ir_led) 

    def shutdown(self):
        # shutdown
        self.i2c.write_byte_data(self.addr, self.mode_cfg, 0b10000000 )

if __name__ == '__main__':
    m = MAX30102()
    for i in range(10):
        time.sleep(1)
        for d in m.read_fifo():
            print(d)

    m.shutdown()

心拍数を計算してみよう

まずはセンサからどんなデータが取得できるかを見てみます。１０秒間計測したデータは図のようになりました。緩やかに数値が上昇し、急速に減少するような山が周期的に見られます。全体的な大きな変動も見られますが、これは装着具合でも変化します。

脈を捉えるためには、山の数がいくつあるかをカウントできれば良さそうです。今は、スマホのセンサを活用した心拍数の計測などもできることもあり、PPG信号のピーク検出アルゴリズムに関する論文は近年も複数見られます。今回は、リアルタイムに処理でき、ノイズに左右されにくい方法として、自己相関関数を用いた方法¹を参考にしました。なお、Raspberry Pi Zeroであれば十分計算できますが、計算能力の低いマイコンではより計算コストの低いアルゴリズムを検討するのが良さそうです。

今回採用したアルゴリズムについて説明します。まず、T秒周期でデータを取得し、自己相関関数を計算します。ここでの自己相関関数は、以下の定義に従います。

f_k = \frac{\sum^{T}_{t=k+1}(y_t-\bar{y})(y_{t-k}-\bar{y})}{ \sum^{T}_{t=1}(y_t-\bar{y})^2}

時系列データ$y$に対し、$\bar{y}$は平均値、$k$は時間のラグを表します。

これにより、脈動したときの周期的なデータの振る舞いが、自己相関関数のピークとして現れます。ここから最初のピーク位置を抽出することで、１回脈動する周期に対応したサンプル数がわかります。図の例では18に最初のピークが現れました。サンプリングレートから、1秒間に25サンプルのデータが生成されるので、25を１周期分のサンプル数18で割り60を掛けた値が、1分間の拍動回数（bpm）になります。

自己相関関数の計算や、ピーク位置の計算は以下のプログラムのとおりです。脈拍数の下限を考慮し２秒周期でデータを取得しようと思います。その場合、平均的には50サンプル程度のデータが取得できるはずなので、自己相関関数で計算するラグの上限値をデフォルトで50としました。ただし、タイミングによってサンプル数が前後するので、50もしくは実際得られたサンプル数のどちらか小さい方をラグの上限値としました。また、ピーク位置を検出する際に、ノイズの影響で出来る小さな山もピークとして検出されるため、orderのオプションを5と設定しました。これは、検出されたピークが、ピーク位置から前後5点の値と比較しても一番大きければピークとみなすという条件です。

import numpy as np

def calc_acf(data, nlag=50):
    y = np.array(data) - np.mean(data)
    nlag = min(len(y),nlag)
    acf = np.array([np.sum(y[lag:]*y[:len(y)-lag]) for lag in range(nlag)])/np.sum(y**2)
    return acf

def find_peak(data,order=5):
    y = np.array(data)
    peak = np.full(len(y),True)
    for o in range(1,order+1):
        diff = y[o:]-y[:len(y)-o]>0
        margin = np.full(o,False)
        peak = np.logical_and(peak,np.logical_and(np.append(margin,diff),np.append(~diff,margin)))
    return np.where(peak)[0]
# dataは２秒間分のデータ
acf = calc_acf(data, fft=False, nlag=50) 
peak = find_peak (acf, order=5) 
bpm = 25/peak[0]*60

消費カロリーを計算しよう

心拍数が取得できたので、これを使って消費カロリーを計算したいと思います。消費カロリーの計算で広く使われる指標が、METs(metabolic equivalents)という単位で、安静時を１としたとき何倍のエネルギーを消費するかの運動強度を表すものです。運動の種類ごとにMETs値は異なり、運動時間と体重をかけ合わせて消費カロリーを計算します。

運動強度は安静時、運動時、最大時のそれぞれの酸素摂取量から算出でき、酸素摂取量は心拍数から算出されることが知られています。そのため、計算には安静時の心拍数や最大心拍数が必要です。最大心拍数は年齢からおおよそ決められますが、安静時の心拍数は朝起きた瞬間の心拍などを記録しておく必要があります。

その他、心拍数と消費エネルギーに関する関係を直接調べた論文もいくつか見つかります。その中で、性別、体重、年齢、心拍数のデータから精度の良い回帰モデルを作成した論文²を参考に計算することにしました。消費エネルギーの計算式は以下の通りです。

$$
{\rm EE} = {\rm Gender} × (-55.0969 + 0.6309 × {\rm HR} + 0.1988 × {\rm Weight } + 0.2017 × {\rm Age}) \
+ (1 - {\rm Gender}) × (-20.4022 + 0.4472 × {\rm HR} - 0.1263 × {\rm Weight} + 0.074 × {\rm Age})
$$

各項目の単位は以下の通りです。

EE(energy expenditure) (kJ/min)
Gender 1 (男性) or 0（女性）
HR(heart rate) beats/min (bpm)
Weight (kg)
Age (years)

なお、1 cal=4.18 Jなので消費カロリー（kcal）を求めるには、EEの値を4.18で割ります。今回はプログラムに直接、性別、体重、年齢を記載するようにしますが、タクトスイッチなどを駆使してデバイス側で設定できるようにすると便利そうですが今後の課題です。上記の計算式も含めた最終的なプログラムは、次の記事に記載します。

ちなみに、これらの消費エネルギーは運動習慣などによっても異なるという研究結果も報告されていますが、ひとまずフィットネスバイクに表示されていたものよりは、尤もらしい目安として活用できると期待してます。

通信方式を考える

今回使用したRaspberry Pi Zero Wと無線通信マイコンとは、WiFiもしくはBluetooth通信が可能です。WiFiは通信距離が長く通信速度も早い特徴がありますが、消費電力が大きい欠点があります。子どもと一緒に遊ぶことを考えると、フィットネスバイクを漕ぐ位置から目の届く範囲でプラレールを動かすので、短距離でも消費電力の少ないBluetooth通信が良さそうです。

Bluetoothには様々なバージョンがありますが、今回はv4.1が対応しているようです。バージョン４以降は、Bluetooth ClassicとBluetooth Low Energy (BLE)の双方が利用できます。前者は大きなデータをやり取りできるのが特徴で、ワイヤレスヘッドフォンなどで広く使われています。後者は、大きなデータをやり取りしない代わりに、接続数が多く、通信距離も長く、消費電力が少なくて済む利点を持っています。プログラム面では、ペアリングさえしてしまえば簡単にシリアル通信できるClassicは容易にリッチな機能を作りやすいのですが、消費電力の観点でどこまで差が出るか検討してみました。

ClassicとBLEの両方の方法で、データを送受信し続けるテストを実施しました。その際に、プラレールのマイコン側の消費電力を、USBの消費電力チェッカーを使って計測したところ次の結果になりました。

Classic方式：10分間で10mA消費
BLE方式：10分経っても0mAのまま変化なし

今回使っている電池の容量が750 mAhなので、Classicでは通信だけでも１時間くらいしか持たないことになります。長く遊ぶにはBLE通信を選択するのが良さそうです。

BLE通信では、セントラル、ペリフェラルという２つの役割に分けられます。なお、Arduino IDEでサンプルコードを利用する際には、サーバー、クライアントと表記されているので混乱ポイントですが、それぞれ以下のような特徴があります。

セントラル：ネットワーク上の親局。データを読み取る側であるためクライアントとも呼ばれる。
ペリフェラル：セントラルに対して情報を提供する子局。データを送信する側のためサーバーとも呼ばれる。

今回は、Raspberr Pi側は回転速度を送るのでペリフェラル、プラレール側はそれを受け取るのでセントラルとして設定しています。大まかな処理の流れは図の通りです。

まず、ペリフェラルがセントラルに見つけたもらえるようにアドバタイジングします
セントラルがペリフェラルをスキャンし、指定するService UUIDを見つけます。
対象の機器を見つけたらペリフェラルに接続要求を行います。
ペリフェラルがセントラルから接続要求を受け接続を完了します。
ペリフェラルがセントラルに向け定期的にデータ送信を行います。今回の通信ではデータが多少欠損しても問題ないので受信完了の確認は不要としました。そこで、データ部分であるキャラクタリスティック（Characteristic UUIDを指定）の属性をNotifyとしました。
終了時には接続を切断します。

おわりに

次回は、いよいよこれまで作成したIoT電池とIoTサイクルメータを電波で繋ぐ実装をしていきます。

電子工作に興味を持ち始め専門ショップに足を運んでみると、多種多様なセンサがあってとてもワクワクします。しかも、一つ一つはお手頃価格なので色々試してみたくなり、つい勢いで買ってしまいますが、簡単なのは買うところまでです。使い方が書かれているデータシートは入ってますが、知識が無いと読み解くのに苦戦を強いられ、プログラム書くよりも使い方を調べる時間の方が長かった気もします（笑）

パルスオキシメータは小さい基板ながらデータシートは割と文量があり、２つのLEDのタイミングの制御など電子工学的な内容や、LEDの波長やフォトダイオードの量子効率など光学的な内容なども関係してくるため、読み解くのは非常に骨が折れました。もちろん、こんな面倒なことをしなくてもライブラリを提供してくださってる有志の皆さんがいらっしゃるので、先人の力を借りればプログラムを実行するだけで難なく心拍数やSpO2が計算できてしまいます。

一方で、「信号強度が弱くてデータが読み取りにくい」「時間が経つに連れ誤差が増えてくる」など取得データに異常があった場合どうなるでしょう。データ処理側でうまく対処することを真っ先に思いつくかもしれませんが、もしデータシートを理解していたら「信号強度を上げるためLEDの電流値を上げる」「熱による誤差を補正するためデバイス温度を取得する」などデータ取得側の対策で精度がグッと上がる可能性が見えてくるのではないでしょうか。

データ分析においても、データを分析する側と取得する側の双方で、原理・現象についての相互理解が深まれば、より高い効果を得られるアイディアが生まれてくる場面も多いのではないかと思います。センサモジュールのデータシート理解は、データサイエンスにおけるデータ理解の第一歩にも繋がるかもしれないので、少しずつ読み解きながら新しい発見を楽しんでみるのも良いかもしれません。

S. Das, S. Pal and M. Mitra, "Real time heart rate detection from PPG signal in noisy environment," 2016 International Conference on Intelligent Control Power and Instrumentation (ICICPI), pp. 70-73 (2016), doi: 10.1109/ICICPI.2016.7859676. ↩
LR. Keytel, JH. Goedecke, TD. Noakes, H. Hiiloskorpi, R. Laukkanen, L. van der Merwe and EV. Lambert. "Prediction of energy expenditure from heart rate monitoring during submaximal exercise", J Sports Sci. 23(3):289-97 (2005), doi: 10.1080/02640410470001730089. PMID: 15966347. ↩

IoTプラレールを作ってみました（その5：プラレールとフィットネスバイクを電波でつなごう）