オプティカルフローの可視化

細胞は変形しながら移動することができます。これを細胞遊走といいます。今回は細胞の移動を可視化、定量化するためのチュートリアルです。

まず、解析の対象とする動画として、Nature誌のLeader cells regulate collective cell migration via Rac activation in the downstream signaling of integrin β1 and PI3Kという論文に含まれるものをダウンロードします。Supplementally information Videos Movie 1を使うことにします。

オプティカルフロー

オプティカルフローとは、物体の動きをベクトルで表すものです。OpenCVではいくつかのアルゴリズムが実装されています。今回は、Farnebackのアルゴリズムを実装した、calcOpticalFlowFarneback関数を使います。難しいオプティカルフローの部分はこの関数が1行で行ってくれるため、プログラムに書くのは、この関数が計算してくれたベクトルを描画するなどの周辺的なことです。オプティカルフローの計算は、2つの画像を元に行います。すなわち、画像を細かい区画に分割し、それぞれで局所的なパターンマッチのようなものを行い、そのマッチ先へのベクトルを計算します。

calcOpticalFlowFarneback(HIS_source, source, flow, 0.8, 10, 15, 3, 5, 1.1, 0);

としていますが、それぞれの引数の意味は以下の通りです。ここで、Mat型の flow はオプティカルフローの計算に一時的に必要なものと思ってください。

cv::calcOpticalFlowFarneback

void calcOpticalFlowFarneback(const Mat& prevImg, const Mat& nextImg, Mat& flow, double pyrScale, int levels, int winsize, int iterations, int polyN, double polySigma, int flags)

Gunnar Farneback のアルゴリズムを用いて,密なオプティカルフローを求めます.

prevImg – 8ビット,シングルチャンネルの1番目の入力画像.
nextImg – prevImg と同じサイズ,同じ型の2番目の入力画像.
flow – 計算されたフロー画像.サイズは prevImg と同じで,型は CV_32FC2 .
pyrScale – それぞれの画像に対する画像ピラミッドを作るためのスケール (<1) を指定します. pyrScale=0.5 は,隣り合う各層において,各層が前の層の半分のサイズになっている古典的画像ピラミッドを意味します.
levels – 最初の画像を含む,画像ピラミッドの層の数. levels=1 は,追加の層が作成されず,元画像だけが利用されることを意味します.
winsize – 平均化窓サイズ.この値が大きくなると,画像のノイズに対するアルゴリズムの頑健さが増し,高速なモーションを検出できる場合が多くなります.しかし,ボケたモーションフィールドを生成することになります.
iterations – 画像ピラミッドの各レベルにおける,アルゴリズムの反復数.
polyN – 各ピクセルにおける多項式展開を求めるために利用される,ピクセル近傍領域のサイズ.この値が大きくなると,画像はより滑らかなサーフェイスで近似され,アルゴリズムの頑健さは増します,しかし同時に,ボケたモーションフィールドも増加します.一般的には, polyN =5 あるいは 7 となります.
polySigma – 多項式展開の基底として利用される導関数を滑らかにするための,ガウス分布の標準偏差. polyN=5 ならばpolySigma=1.1 , polyN=7 ならば polySigma=1.5 が適当な値です.
flags –
処理フラグ.以下の値の組み合わせ:

OPTFLOW_USE_INITIAL_FLOW 入力 flow をフローの初期推定値として利用します.
OPTFLOW_FARNEBACK_GAUSSIAN オプティカルフロー推定のために,同サイズのボックスフィルタの代わりに winsize*winsize サイズのガウシアンフィルタを利用します.通常,このオプションによって処理速度は低下しますが,ボックスフィルタを利用する場合よりも正確なフローを求めることができます(また,同程度の頑健性を得るためには,ガウシアン窓の winsize をより大きくしなければいけません).
[ opencv 2.2 documentation より]

細胞遊走の動画を入力に用いると次のように各局所区画の移動方向が線分で表示されます。
opti-s

opti-o
#include "stdlib.h"
#include "opencv/cv.h"
#include "opencv/highgui.h"
#include "string.h"
#include "stdio.h"
#include "iostream"
#include "math.h"

using namespace std;
using namespace cv;

int main(){

    VideoCapture cap("srep07656-s2.avi");
    
    int Width=cap.get(CV_CAP_PROP_FRAME_WIDTH);
    int Height=cap.get(CV_CAP_PROP_FRAME_HEIGHT);
    int max_frame=cap.get(CV_CAP_PROP_FRAME_COUNT); //フレーム数
    
    
    Mat source(Height,Width,CV_8UC1);
    Mat HIS_source(Height,Width,CV_8UC1);
    
    for(int frame=0; frame<max_frame; frame++){
        cap>>source;
        Mat disp=source.clone();
        cvtColor(source, source, CV_BGR2GRAY);
        
        if(frame>0){
            
            vector<cv::Point2f> prev_pts;
            vector<cv::Point2f> next_pts;
            
            Size flowSize(100,100); //ベクトルの数
            Point2f center = cv::Point(source.cols/2., source.rows/2.);
            for(int i=0; i<flowSize.width; ++i) {
                for(int j=0; j<flowSize.width; ++j) {
                    Point2f p(i*float(source.cols)/(flowSize.width-1),
                                  j*float(source.rows)/(flowSize.height-1));
                    prev_pts.push_back((p-center)*0.95f+center);
                }
            }
            
            Mat flow;
            vector<float> error;
            
            calcOpticalFlowFarneback(HIS_source, source, flow, 0.8, 10, 15, 3, 5, 1.1, 0);
            
            // オプティカルフローの表示
            std::vector<cv::Point2f>::const_iterator p = prev_pts.begin();
            for(; p!=prev_pts.end(); ++p) {
                const cv::Point2f& fxy = flow.at<cv::Point2f>(p->y, p->x);
                cv::line(disp, *p, *p+fxy*8, cv::Scalar(0),1);
            }
            
            HIS_source=source.clone();
            
            imshow("vector", disp);
            imshow("source", HIS_source);
            
            int c = waitKey(1);
            if(c==27)return 0;
        }
        
        cout<<frame<<endl;
        frame+=1;
    }
    waitKey(0);
    return 0;
}