最も近い値を持つデータを別のデータセットからマッチングする方法について考える話

まずデータセットを2つ

data Q1;
call streaminit(777);
do ID=1 to 40;
X = round(rand('uniform')*1000,0.1);
output;
end;
run ;

data Q2 ;
call streaminit(888);
do ID2=1 to 40;
Y = round(rand('uniform')*1000,0.1);
output;
end;
run ;

それぞれのデータセットについて XとYの値はデータセット内でたまたま一意であることを確認しています。

で、やりたいことは

Q1のXの値に最も近いYの値を持つQ2のデータをQ1の各オブザベーションに結合するということです。

最も近いというのはX-Yの絶対値が最小になるデータです。

求めるべき結果は

さて、どうするかという話で、まず思いつくのはsqlやデータステップで

今回の場合、40obs×40obs=1600obsの総当たりの直積をつくってIDとabs(X-Y）でソートして

first.IDをとるような処理ですね。

それもありなんですが、例えばデータが1000や10000、さらにもっと増えた場合、処理時間がえらいことになっていくので今回は別の方法を検討してみます。

まずは、

data _Q2;

   set Q2 end=eof;

   CLUSTER=_n_;

   if eof then call symputx('nobs',_N_);

run;

proc fastclus data=Q1 out=_Q1

              seed=_Q2(rename=(Y=X))  maxclusters=&nobs

              noprint maxiter=0 ;

     var X;

run;

proc sort data=_Q1;

by cluster;

run;

proc sort data=_Q2;

by cluster;

run;

data A2;

informat ID X ID2 Y;

     merge _Q1(in=ina) _Q2;

     by cluster;

if ina;

keep ID X ID2 Y distance;

run;

で先ほどの結果になります。

何をしているかというと、クラスター分析ができるfastclusプロシジャを使って、Q2のYをシード値にして、最大クラスタ数をQ2の総オブザベーション数にします。

あとはclusterでマージすれば終わりです。

続きましては

data A1;

  if _N_ = 1 then do ;

if 0 then set Q2;

  declare hash h1 (dataset:'Q2', ordered:'A') ;

  declare hiter hi ('h1');

  h1.definekey('Y');

h1.definedata('Y','ID2');

  h1.definedone () ;

end ;

  set Q1;

rc = h1.find(key:X);

  if rc ^=0 then do;

  h1.add(key:X,data:X,data:ID);

  if hi.setcur(key:X) = 0 and hi.prev() = 0 then do;PID = ID2;PY = Y ;end;

  if hi.setcur(key:X) = 0 and hi.next() = 0 then do;NID = ID2;NY = Y ;end;

if PY=. and NY=. then put ID=;

if PY =. then do;Y=NY;ID2=NID;end;

else if NY =. then do;Y=PY;ID2=PID;end;

else if abs(X-PY)<=abs(X-NY) then do;Y=PY;ID2=PID;end;

else if abs(X-NY)<abs(X-PY)  then do;Y=NY;ID2=NID;end;

  h1.remove(key:X);

end;

distance = abs(X-Y);

run ;

で、適当に

data A1;

informat ID X ID2 Y;

set A1;

keep ID X ID2 Y distance;

run;

成型すれば先ほどの結果と同じです。

ハッシュオブジェクトなわけなんですが、何が面白いかって、総当たり比較を行っていないってことなんですね。ハッシュオブジェクトにQ2を入れて、データステップでXを読み込んで、まずfindメソッドで差が0のものがないかをみます。

それがなかった場合、X自身をハッシュオブジェクトに入れます。ハッシュオブジェクトは定義時に値で昇順に格納されるようにしているため、X自身が格納された場所の前のデータと後のデータを取得することで、最近値候補の2つをゲットできます。あとはその2つとの差を比べて、どっちを選択するかだけです。

setcurメソッドは指定したキー値の値で、ハッシュ反復子オブジェクト内のポインタをそこに持っていきます。prevメソッドとnextメソッドはポインタの前後のデータを取得するメソッドです。

処理が終了したらハッシュにいれたXはremoveで取り除いておきます。

個人的に勉強中で、謎なのがremoveメソッドは通常、ハッシュ反復子オブジェクトで再定義されているハッシュオブジェクトを対象にして、かつ以前にハッシュ反復子オブジェクト対象のメソッドが実行されていた場合、ロックがかかってエラーになるはずなんですか、

今回のようにif 条件式の中で判定のために使うと、メソッドの結果を得つつ、ロックを回避できるんですよね。

もし、このことについて知っている方がいらっしゃれば情報求むです。
※後日追記：コメントでamatsuさんに突っ込んで貰いました。僕の勘違いです。条件式うんぬんではなくハッシュ反復子のポインタの位置があるキーをremoveできないのです。アホデシタ。

今はオブザベーション数40ですが、これを1000や10000にしてみると（今回はhashにmultikeyつけてないので値は一意になるように）、ハッシュオブジェクトの処理時間が相当短いことが確認できます(SAS雲丹だと厳しいけど）。さらに10万や100万規模になってもhashexpでハッシュオブジェクトのサイズを調整すれば充分優秀です。