コロナ肺炎とパルスオキシメーター
新型コロナウイルスの重症度の判定で血中酸素飽和度を簡単に測定できるでバルスオキシメーターが品薄になってしまった。パルスオキシメーターを指に挟むことで脈拍(脈波・心拍数)と酸素飽和度(サチュレーション)を簡単に測定できる。歯科治療においても、この2項目をモニターすることで、患者さんの異常を早く検知できるのでとても心強いし、価格も安価(だった)なので複数持っている先生も多いと思う。
体温計が一気に市場から消えたのと同じく、新型コロナウイルスによる肺炎の重症度の目安として「血中酸素飽和度(SpO2)が常に95%以下」との報道がなされると同時に店から消えた。中国発送のものはまだ買えるみたいだけど、到着まで早くても1か月以上はかかると思う。
昔にラズパイのセンサーを買い漁ってた時にSpO2を測定できるセンサーのMAX30102を2個買っていたことを思い出した。この時は、MAX30102をラズパイで使うライブラリを見つけられなくて早々に興味を無くし、いつか使うであろうボックスに投げ込んでおいた。
このパルスオキシメーターが足りない現在、ArduinoやESP8266で実用化できるレベルまで使えるようにすれば、もしコロナに感染した時に常にモニターして、心拍、SpO2に異常があれば家族に自動的にメールなどすれば最悪の事態は避けられるかもしれない。
MAX30102(心拍数・SPO2センサー)をESP8266で使う
早速、ESP8266、MAX30102で下のサイトを参考にすすめてみた。
I2C接続でMAX30102を繋げる。Arduinoのデフォルト(A4[SDA]・A5[SCL])で接続することになるけど、ESP8266で利用するにはGPIO0とGPIO2を使って行う。基本構成は下図で行こうとおもう。ESP-01な感じで書いてあるけど( ´∀｀ )。。。技適通ったESP8266に適宜読み替えてください。
Wire.begin(0,2);
をI2Cセットアップ時にを忘れずに！！。
ちなみにMAX30102はしたのを購入済み(去年の10月w)
で、参考サイトを読んでいると、ライブラリも導入する必要がある。ライブラリは下
Arduino Studioからもインストール可能。サンプルスケッチExample8_SPO2を実行してみる
ん・・・。HR=166、SPO2=76-89とか・・・。
このバイタルサインは・・・死んでますよ？ｗｗ
どんな処理しているのかコードを追ってみる。。。どうやらspo2_algorithm.cppが心臓部らしい。
読んでみると・・・
なにしてるかよくわかんねぇぇよｗｗ
「これは・・・終わったかな？」とおもったけど、下のテストスケッチで生データをシリアルプロッタでみると、ちゃんとデータは取れているようだ・・・。
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 |
#include <Wire.h> #include "MAX30105.h" MAX30105 particleSensor; void setup() { Wire.begin(0,2); // I2Cのピン設定 Serial.begin(115200); Serial.println(""); Serial.println("Initializing..."); // Initialize sensor while(!particleSensor.begin(Wire, I2C_SPEED_FAST)){ Serial.print("."); } Serial.println("OK!"); //Setup to sense a nice looking saw tooth on the plotter byte ledBrightness = 0x1F; //Options: 0=Off to 255=50mA byte sampleAverage = 8; //Options: 1, 2, 4, 8, 16, 32 byte ledMode = 2; //Options: 1 = Red only, 2 = Red + IR, 3 = Red + IR + Green int sampleRate = 400; //Options: 50, 100, 200, 400, 800, 1000, 1600, 3200 int pulseWidth = 411; //Options: 69, 118, 215, 411 int adcRange = 4096; //Options: 2048, 4096, 8192, 16384 particleSensor.setup(ledBrightness, sampleAverage, ledMode, sampleRate, pulseWidth, adcRange); //Configure sensor with these settings Serial.println("ir,red"); } void loop() { uint32_t ir_v = particleSensor.getIR(); uint32_t red_v = particleSensor.getRed(); Serial.print(ir_v); Serial.print(","); Serial.println(red_v); } |
赤外線(IR)パルスと赤色(Red)パルスを分析することでSPO2がわかるらしいので、できるかどうかも含めて資料をネットの海から集める。
https://arxiv.org/ftp/arxiv/papers/1907/1907.11989.pdf
(抜粋)上のグラフ図にすべてが集約されているようだ。つまりはRed_AC、Red_DC、IR_AC、IR_DCがわかればSPO2が出るってことね。あとは
$$R=\frac{Red_{AC} \div Red_{DC}}{IR_{AC} \div IR_{DC}}$$
でRを求めて、RからSPO2を求める式はspo2_algorithm.cppにあったのでこれを準用する。
$$SPO_{2}[\%]=-45.060R^2+30.354R+94.845$$
自分はSPO2は通常、安静状態で97-99程度(手持ちのパルスオキシメーターで確認済)なので、
$$98=-45.060R^2+30.354R+94.845$$
上の二次方程式を解くと 大体R = 0.1284 , 0.5452となる。この付近の数値になればいいと予想がつく。
その上でグラフの目盛りを大体で読んで手計算で求めたら、R=1.782とかになった。何度も何度も検算しても1以下にならない
もうさぁ～このセンサー壊れてない？
心が折れそうになったその時、1.782って・・・なんか・・・
$$\frac{1}{1.782}=0.561$$
逆数にするとちょうどよくね？息ごらえ(SPO2を意図的に下げる)して計算しても93-94%に落ちているようだ。
これは・・・IRとRedが逆なんでは？
コード上でが逆になってる様子はない。嫌な予感・・・ledMode = 1(赤色LEDのみ)にしてセンサーをみると、赤色LEDは光っていない。指をかざすとgetRed()で変化する数値が取得で来ている・・・・これは・・・もしや・・・IR(赤外線)のほうでは？
赤外線カメラで確認してみると
おEEEE！！！やっぱり逆じゃねーか！
この個体が逆なのか、メーカ正規品(SparkFun)は正常なのか、わからないけど、これはとりあえずIRとRedは逆なんですねw。こんなので2日も悩んだよ・・・マジで・・・
組み立てとスケッチ実装
気を取り直して、3Dプリンタでマウンタを作り、マジックテープのバンドをつけて安定して固定できるようにした。
センサー表面には直接金属端子が触れないよう、セロテープを表面に張り付けた。
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 |
#include <Wire.h> #include "MAX30105.h" #include <movingAvg.h> MAX30105 particleSensor; //------ 定数定義 const uint32_t TH_FIN = 7000; //指が置いてあるかどうかの閾値 const int32_t TH_AMOUNT = 300; //パルスの起点となる検出変化量 const int32_t MIN_INIT = 9999999; //最小値の初期値 const int32_t MAX_INIT = 0; //最大値の初期値 //表示する心拍数とSpO2の範囲(用途により適宜変更) const uint32_t DISP_MIN_HR = 30; const uint32_t DISP_MAX_HR = 180; const uint32_t DISP_MIN_SPO2 = 70; const uint32_t DISP_MAX_SPO2 = 100; //パルス検出保持変数 long l_time = millis(); //最後のパルス検出時間の保持 int32_t before_ir_v = 0; //最後のir_v値保持 int32_t b_diff = 0; // 最後の差の保持 long pulse_interval = -1; // パルス間隔時間 //1パルス中の生データの最大・最小値の保持 int32_t min_ir_v = MIN_INIT ,max_ir_v = MAX_INIT; int32_t min_red_v= MIN_INIT ,max_red_v = MAX_INIT; //移動平均値(IR_DC、RED_DC) 30(20～50)サンプル movingAvg avgIr_v(30); movingAvg avgRed_v(30); //移動平均値(心拍数は直近3ビート、SPO2は5ビートの平均) movingAvg avgHR(3); movingAvg avgSPO2(5); void setup() { Wire.begin(0,2); // I2Cのピン設定 Serial.begin(115200); Serial.println(""); Serial.println("Initializing..."); // センサー初期化 while(!particleSensor.begin(Wire, I2C_SPEED_FAST)){ Serial.print("."); } Serial.println("OK!"); byte ledBrightness = 0x1F; //Options: 0=Off to 255=50mA byte diffmpleAverage = 8; //Options: 1, 2, 4, 8, 16, 32 byte ledMode = 2; //Options: 1 = Red only, 2 = Red + IR, 3 = Red + IR + Green int diffmpleRate = 400; //Options: 50, 100, 200, 400, 800, 1000, 1600, 3200 int pulse_intervalWidth = 411; //Options: 69, 118, 215, 411 int adcRange = 4096; //Options: 2048, 4096, 8192, 16384 particleSensor.setup(ledBrightness, diffmpleAverage, ledMode, diffmpleRate, pulse_intervalWidth, adcRange); //Configure sensor with these settings //初期値として一回読んでおく(ir_vをパルス計測に用いる) before_ir_v = particleSensor.getRed(); //パルス最終取得時間の初期化 l_time = millis(); //移動平均ライブラリの初期化 avgIr_v.begin(); avgRed_v.begin(); avgHR.begin(); avgSPO2.begin(); } void loop() { //逆になっている。センサーが逆を返しているので各自チェック uint32_t red_v = particleSensor.getIR(); uint32_t ir_v = particleSensor.getRed(); //指を置いているか if(red_v< TH_FIN || ir_v < TH_FIN) return; //移動平均値(IR_DC RED_DCの算出) double ir_v_dc = avgIr_v.reading(ir_v); double red_v_dc = avgRed_v.reading(red_v); //IRとREDのACを求める為に最大値・最初値の更新 if(ir_v<min_ir_v) min_ir_v = ir_v; if(ir_v>max_ir_v) max_ir_v = ir_v; if(red_v<min_red_v) min_red_v = red_v; if(red_v>max_red_v) max_red_v = red_v; //パルス検出にはir_vを利用。前回との差分(変化値)を求める int32_t diff = before_ir_v - ir_v; //TH_AMOUNT以上の変化がある場合、パルス起点とする if(b_diff < TH_AMOUNT && diff > TH_AMOUNT){ //1パルスの時間を計算 pulse_interval = millis() - l_time; l_time = millis(); //1パルスの時間より心拍数の計算。avgHRは整数のみなので1000倍して1000で割り少数保持 double hr = (double)avgHR.reading(60000*1000/pulse_interval) / 1000.0; //SPO2の計算 //IR・REDのACを求める(振幅最大-振幅最小) int32_t ir_v_ac = max_ir_v-min_ir_v; int32_t red_v_ac = max_red_v-min_red_v; // R = (AC_RED / DC_RED) / (AC_IR / DC_IR)の計算式より double red_div = double(red_v_ac)/red_v_dc; double ir_div = double(ir_v_ac)/ir_v_dc; double R = red_div / ir_div; // SPO2 = -45.060*R^2 + 30.354*R + 94.845 これはspo2_algorithm.cppにあったのを準用 // 乗除1000は少数保持のため double spo2 = (double)avgSPO2.reading((-45.060*R*R + 30.354*R + 94.845)*1000.0) / 1000.0; //最大値・最小値の初期化 min_ir_v = MIN_INIT; max_ir_v = MAX_INIT; min_red_v = MIN_INIT; max_red_v = MAX_INIT; //心拍数とSPO2の表示(範囲は定義で変更可) if(hr <= DISP_MAX_HR && hr >= DISP_MIN_HR && spo2 <= DISP_MAX_SPO2 && spo2 >= DISP_MIN_SPO2){ Serial.print(hr); Serial.print(","); Serial.println(spo2); } } //パルス検出用値の保持 before_ir_v = ir_v; b_diff = diff; } |
スケッチにコメントいれてあるけど、要点は3つのみ
- ビートの区切りは直前のデータとの差の変位量で判断(IRの方を使っている)
- DCの計算は1ビートの移動平均を使う
- ACの計算は1ビートの最大値と最小値の差で計算
ちなみにビートとDC(移動平均)のグラフを掲載
ちなみに移動平均を簡単に処理してくれる下のライブラリを使わせてもらった。感謝！
ノイズ処理はもーちょいしたほうがよさそうだけど、サンプルスケッチなんかよりはよっぽど安定していると思う。
赤はSPO2、青は心拍数。SPO2が下がっているところは息ごらえ(頑張って息を止める)をしたところ。
酸素が少ない血が指まで運ばれるのに時間がかかるので、実際下がるのは5秒～10秒後なので下がらないからといって止めてると死にますから注意です。
心拍数はギザギザ上下しているように見えけど、これは息を吸う時は心拍数が少し上がり、吐くときに少し下がるという生理現象。さすがに息止めの前後は変化してるけど、呼吸性の変化も捉えることができているのは驚き。
手持ちのパルスオキシメーターとの比較したけど、ほぼ同じ数値をとっていました～。とりあえず、ここまでできればあとはデータベースで焼くなり煮るなりすればいいと思います。95%以下が3分続いたらメールとかね。もう疲れたから今度やります。
最後に一言！センサーを使う前に
赤外線と赤色光の数値を要チェック！！
これ重要！！
[追記2020/05] リアルタイムにグラフ表示もしてみました
[追記2020/05/20]
FBでライブラリのissueに「IRとRedが逆じゃない？」みたいな質問あるねって情報もらいました。
5/10にissue投稿してるんだけど、5/4-5にこの記事だしてるから、もしかしてこのサイトみたのかなぁ(期待)？それはいいとして…とりあえずSparkFun純正使ってないと答えてくれないだろうし…回答を見守ります。
[追記2021/11/05]
ぶらぶらしていたら結構本格的に解説しているサイトがありました。いまちょっと別なことやってるんで今度読み込んでみたいと思います
