Word eXtraction method、いわゆる文法圧縮の一種。開発者はOkanohara氏。20年以上前の化石programですが、面白いので紹介します。
概要
WX法は、Dataを最小構成要素(以下これをUnitと呼ぶ)に分解する算法です。
Unitは例えば、自然言語で言えば、単語などに相当し、DNAの配列情報で言えば遺伝暗号(codon)に相当します。
様々なDataに対してそれぞれのUnitを1つ1つ定義するのは費用がかさむ上、柔軟性がありません。そこで、Unitは、その内部に、さらに小さい構成要素を持たないという定義とします。
WX法はDataから統計情報を用いて、そのような独立した部分列を抽出し、全Dataを上の評価基準で尤も的確な形に分解する算法です。
WXの算法
WX法は次の5つの手順からなります
- 部分列を列挙する
- それぞれの部分列に対し、Unitらしさを表す評価値を与える
- (2)の評価値を用いて、Unit分解を行う
- (3)の結果を用いて評価値を再計算する
- 収束するまで(3)(4)を繰り返す
- 最終的にentropy符号で符号化
性能
圧縮率は大きなfileに対してLHAやgzipを超え得る。古典的なLZWも概ね凌駕します。PPMやBWT系には余裕で負けるでしょう。
圧縮速度は遅く、展開は高速。
実装編
本家はC++ですが、例によってJavaScriptに翻訳。ついでなので改良もしておきます。主な変更点は…
- 圧縮文に刻印する辞書自体も圧縮
- 頻度表も圧縮
- 静的Range Coderだけど符号化途中で頻度表を更新
の3点です。そのせいで若干低速化しています。
const culcPoint=(w,l,fe,le)=>w?0|Math.log(w)*l*1024:1<<31,entropyTable=function(s,S,a,c,d,b,i,p,w){for(;--i||(b=a.charCodeAt(d++),i=a.charCodeAt(8+d));)S[S(b++)]=c++;for(a=i='';c=s.charAt(i++);)if((d=S[c])<92)a+=S(b+d);else if(d<1976){if(p=(w=a.slice(-28)).substring(d=(c=d-92)%28,c=d+c/28+2))for(;w.length<c;)w+=p;a+=w.substring(c,d)}else d<2068?b=12353:d&&(b=12449),a+=c;for(s=[i=0];i<4095;)s[++i]=s[i-1]+a.charCodeAt(i)-64;return s}('@HGFEEDEDCCBA@?@',String.fromCharCode,' ?\u0758[]',0,0);
let N_EXT=16;
function nextMatch(p,P,B){
for(var i=0,a=P[p],b=P[p+1];i<N_EXT&&B[a+i]===B[b+i];)i++;
return i
}
function addWord(a,z,l,B,C,P,D){
let w=z-a,fe=0,le=0;
if(l^1){
if(w<=D.widthBorder)return;
let i=a||1,c;
//entropy検査
for(;i<z;C[B[c+l]+256]++)c=P[i++],C[B[c-1]]++;
//entropyを計算
for(i=256;i;C[i]=C[i+256]=0){
c=C[--i]<<12;
if(c>=w)fe+=entropyTable[c/w-1|0];//-1はentropyTable が 0から4095までだから
c=C[i+256]<<12;
if(c>=w)le+=entropyTable[c/w-1|0];
}
if(fe+le<3000)return
}
D.push({first:a,last:z,length:l,point:culcPoint(w,l,fe,le),number:D.length+1,fe,le,rcount:0,stop:0})//+1なのは1から番号開始のため。0は特別なflag
}
function dicWord(n,D,B,P,O){
for(let i=0,l=D[n].length,p=P[D[n].first];i<l;)
O[O.p++]=B[p+i++]
}
function occupy(n,D,W,P){
for(var{first,last,length,number}=D[n],w=W.length,p,rc=0,i=first,j;i<last;){
p=P[i++];
//空いているか検査.ここは高速化可能
if(p+length>w)continue;
for(j=0;j<length&&!W[p+j];)j++;
if(j<length)continue;//占有されていた
//占有flagを入れる
rc++;W[p]=number;//開始位置
for(j=0;j<length;)W[p+j++]=-number//開始でない場所には負の登録番号をいれておく
}return rc
}
function culcFileSize(D){
let fs=0,ds=0;//file内容の圧縮後の尺, 辞書内容の圧縮後の尺
//圧縮後file尺を予測する
const LOG2=Math.LN2,sumlog2=Math.log(D.wordN)/LOG2;
for(let i=D.length,c;i;)
if(c=D[--i].rcount)
fs+=c*=sumlog2-Math.log(c)/LOG2,
ds+=D[i].length+3;//仮に長さ1byte,確率2byte
return fs/8+ds>>>0
}
function allocate(D,W,P){
let kn=0,wn=0,i=W.length,ds=D.length;
for(;i--;)W[i]=0;
//割り当て
for(;++i<ds;){
if(D[i].stop){D[i].rcount=0;continue}
let c=D[i].rcount=occupy(i,D,W,P);
if(c>0)wn+=c,kn++;
}
D.kindN=kn;D.wordN=wn;
//return culcFileSize(D)
}
function allocate2(D,W,P){
//それぞれのmaxを決める
for(var i=W.length,w=i,l=D.length,n=0,n2i=new Uint32Array(l+1);i;)W[--i]=0;
for(;i<l;){
let{first,last,number,length}=D[i],j=first,p;
for(D[n2i[number]=i++].rcount=0;j<last;)
if(!W[p=P[j++]]&&p+length<=w)W[p]=number
}
//頭から順番に見ていく
for(i=D.wordN=0;i<w;n++){
let d=D[n2i[l=W[i++]]];
d.rcount++;
for(d=d.length;--d;)W[i++]=-l
}
D.wordN=n//;return culcFileSize(D)
}
function extract(B,W,P,z){
var D=[],k=N_EXT,i=k-1,guardPos=z+i,c=B[0],C=new Uint32Array(65536),p;
for(;i;)B[z+--i]=B[i];
for(;i<z;C[c<<8|(c=B[++i])]++)P[i]=i;
for(c=C[i=0];c<z;)c=C[++i]+=c;
//radix sort
for(;k>0;){
for(k-=2,i=z;i;W[--C[B[p+k]<<8|B[p+k+1]]]=p)p=P[--i];
for(k-=2;i<z;P[C[B[p+k]<<8|B[p+k+1]]++]=p)p=W[i++]
}
B[P[z]=guardPos]=B[P[z-1]]+1&255;//最後のぎりぎり
C.fill(0,0,512);D.widthBorder=2;
//単語を抽出. W:前回一致した所の開始位置を保持
for(let pn=i=0;i<z;i++){
c=nextMatch(i,P,B);//次の部分列との一致数を計算
if(c>pn)//前回より多く一致
for(k=c;k>pn;)W[--k]=i;//登録用に今の位置を記録
else if(c<pn)//前回より少なく一致
for(k=c;k<pn;)addWord(W[k],i+1,++k,B,C,P,D);//順次登録していく。登録内容は first<=範囲<last length:長さ
else c||pn||addWord(i,i+1,1,B,C,P,D);//1文字しか出現しないものは特別に登録
pn=c
}
const ds=D.length,cmp=(a,b)=>b.point-a.point;
qsort(D,0,ds-1,"point"),allocate2(D,W,P);
for(k=2;k--;qsort(D,0,ds-1,"point"),allocate2(D,W,P))
for(i=ds;i;c.point=culcPoint(c.rcount,c.length,c.fe,c.le))c=D[--i];//point再計算
for(qsort(D,0,ds-1,"rcount");i<ds&&D[i].rcount;)i++;
D.kindN=i;return D
}
//write counter
function hits(O,H,k,n){
for(var i=0,a=1,b=2,c=k,p=k,F=new Float64Array(33),S=new Float64Array(33),N=n,M=n;p;)F[log2b(H[--p]+1)]++;
for(;i<33&&c>0;c-=f){
var f=F[i++],l=M/(b-1)>>>0,h=c-(c*(b-2)<N);
O.ev(l<c?c-l:0,f>>>0,h>>>0);
N-=f*(b-2);M-=f*(a-1),a*=2,b*=2
}
for(i=33,a=(1<<30)*4,b=a*2,c=N=M=0;i--;b/=2)
S[i]=c,c+=f=F[i],N+=f*(a-1),M+=f*(b-2),a/=2;
for(;n;n-=k,O.ev(l>>>0,k>>>0,h>>>0)){
for(i=0,a=1,b=2,k=H[p++],c=log2b(k+1);i<c;S[i++]--,a*=2,b*=2)
if(f=F[i])O.enc(f>>>0,S[i]>>>0,f+S[i]>>>0);
S[i]>0&&O.enc(0,F[i]>>>0,F[i]+S[i]>>>0);
F[i]--;N-=--a;M-=b-=2;
(l=n-M)<a&&(l=a);
(h=n-N)>b&&(h=b)
}
}
//math.log2
function log2b(a,b,c){
c=(65535<a)<<4,c|=b=(255<(a>>>=c))<<3,c|=b=(15<(a>>=b))<<2;return(c|=b=(3<(a>>=b))<<1)|a>>b>>1
}
function isort(A,a,z,e){
for(var i=a,b,c,d;a<z;A[c]=b)
if((b=A[c=d=++a])[e]>A[i][e])for(;A[c]=A[--c],c>i;);
else for(;b[e]>A[--d][e];)A[c]=A[c=d]
}
function qsort(A,a,z,e){
if(z-a<9)return isort(A,a,z,e);
var b=A[a][e],c,i=A[z][e],j=a,k=z,p=A[a+z>>1][e];
b>i?i>p?p=i:b>p||(p=b):b>p?p=b:p>i&&(p=i);
for(i=a;j<=k;)
if((b=A[j])[e]<p)A[j]=A[k],A[k--]=b;
else{if(b[e]>p)c=A[i],A[i++]=b,A[j]=c;j++}
a<--i&&qsort(A,a,i,e);++k<z&&qsort(A,k,z,e)
}
// compressor
function WXe(A){
function enc(c,f,t){
R=R/t>>>0;L+=c*R;
for(R*=f;R<16777216;L=L<<8>>>0,R*=256)
if((c=0xffffffff<L)||255>L>>>24)
for(O[o++]=255&c--+M,M=L>>>24;N;N--)O[o++]=255&c;
else++N
}
function pb(b,c){
for(L+=(R>>>=1)*b;R<16777216;L=L<<8>>>0,R*=256)
if((c=0xffffffff<L)||255>L>>>24)
for(O[o++]=255&c--+M,M=L>>>24;N;N--)O[o++]=255&c;
else++N
}
function eb(i,b,P,s=5){
var p=P[i],c=(R>>>12)*p;
b?R=c:(R-=c,L+=c);
for(P[i]-=(p-(b<<12)>>s)+b;R<16777216;L=L<<8>>>0,R*=256)
if((c=0xffffffff<L)||255>L>>>24)
for(O[o++]=255&c--+M,M=L>>>24;N;N--)O[o++]=255&c;
else++N
}
function ev(n,v,m,o){
for(;m-n>65535;v>o?n=o+1:m=o)o=n+(m-n>>1),pb(v>o);
n^m&&enc(v-n,1,m-n+1)
}
var z=A.length,B=new Uint8Array(z+N_EXT),P=new Int32Array(z+1),W=new Int32Array(z>N_EXT?z:N_EXT),
o=-1,O=[],L=0,R=o>>>0,M,N=0;
B.set(A);O.eb=eb;O.ev=ev;O.enc=enc;
let D=z<N_EXT||extract(B,W,P,z),i=D.length,ds=i,wn=0|D.wordN,kn=0|D.kindN,s=0,cs=0,E=kn>256?kn+1:257,F=new Uint32Array(E),C=new Uint32Array(E+1),I=new Uint32Array(ds+1);
//write sizes
for(let c=64,n=kn;n-=O[++o]=n&c-1;c=256)n/=c,O[0]+=64;
for(let c=64,n=wn,p=o+1;n-=M=O[++o]=n&c-1;c=256)n/=c,O[p]+=64;
if(z<N_EXT){for(let c of A)O[++o]=c;return O}
for(;i;)I[D[--i].number]=i;
//scan words and write its length
for(let G=Uint16Array.from(Array(60),a=>2048),a=0;i<kn;){
let l=D[i].length,b=l;
//length is gamma coded
for(s+=l;b>>=1;)eb(a++,0,G,5);
for(eb(a,1,G,5);a;)eb(++b+29,l>>--a&1,G,6);
for(b=P[D[i++].first];l;)++F[B[--l+b]]
}
hits(O,F,256,s);
for(i=0;i<s;)C[cs+1]=i+=F[cs++];
//write words
for(E=F.slice(i=0);s;)
for(let l=D[i].length,p=P[D[i++].first],a=0,b;a<l;){
enc(C[b=B[a+++p]],F[b],s);
if(!--E[b])for(b=s=0;b<cs;)C[b+1]=s+=F[b]=E[b++]
}
for(i=0;i<kn;)C[i+1]=s+=P[kn-i]=F[i]=E[i]=D[i++].rcount;
for(P[i=0]=F[kn-1];i<kn;)P[++i]=P[i+1]-P[i];
//write model
for(let G=Uint16Array.from(Array(60),a=>2048),a=i=0;i<kn;){
let l=P[i++],b=l;
eb(29,l>0,G,4);
if(l){
for(;b>>=1;)eb(a++,0,G,5);
for(eb(a,1,G,5);a;)eb(++b+29,l>>--a&1,G,6)
}
}
//write CFG
for(i=0;s;)if((D=W[i++])>0){
enc(C[D=I[D]],F[D],s);
if(!--E[D])for(D=s=0;D<kn;)C[D+1]=s+=F[D]=E[D++]
}
for(;s++<5;L=L<<8>>>0)
if((c=0xffffffff<L)||255>L>>>24)
for(O[o++]=255&c--+M,M=L>>>24;N;N--)O[o++]=255&c;
else++N;
return O
}
// decompressor
function WXd(A){
function cum(t){
for(;R<16777216;R*=256)L=(L<<8|A[a++])>>>0;
R=R/t>>>0;return L/R>>>0
}
function gb(b){
for(;R<16777216;R*=256)L=(L<<8|A[a++])>>>0;
R>>>=1,b=L-R>>>31,L-=R&--b;
return-b
}
function db(i,P,s=5){
for(;R<16777216;R*=256)L=(L<<8|A[a++])>>>0;
let p=P[i],c=(R>>>12)*p,b=1;
L<c?R=c:(R-=c,L-=c,b=0);
P[i]-=(p-(b<<12)>>s)+b;
return b
}
function dv(n,m,o){
for(;m-n>65535;gb()?n=o+1:m=o)o=n+(m-n>>1);
if(n^m)n+=o=cum(m-n+1),L-=R*o;
return n
}
function gets(H,k,n){
for(var i=0,a=1,b=2,c=k,p=0,F=new Float64Array(33),S=new Float64Array(33),N=n,M=n;i<33&&c>0;c-=f){
var l=M/(b-1)>>>0,h=c-(c*(b-2)<N),f=F[i++]=dv(l<c?c-l:0,h>>>0);
N-=f*(b-2);M-=f*(a-1),a*=2,b*=2
}
for(i=33,a=(1<<30)*4,b=a*2,c=N=M=0;i--;b/=2)
S[i]=c,c+=f=F[i],N+=f*(a-1),M+=f*(b-2),a/=2;
for(;n;n-=H[p++]=dv(l>>>0,h>>>0)){
for(i=0,a=1,b=2;S[i]>0;S[i++]--,a*=2,b*=2)if(f=F[i]){
if(cum(c=f+S[i]>>>0)<f){R*=f;break}
L-=R*f;R*=S[i]
}
F[i]--;N-=--a;M-=b-=2;
(l=n-M)<a&&(l=a);
(h=n-N)>b&&(h=b)
}
}
var a=0,i=A[a++],s=64,kn=i&63,wn,L,R=-1>>>0,o=0,O=[];
for(i>>=6;i--;s*=256)kn+=A[a++]*s;//read dict size
i=A[a++],s=64,wn=i&63;
for(i>>=6;i--;s*=256)wn+=A[a++]*s;//read CFG size
if(!kn&&!wn){for(s=A.length;a<s;)O[o++]=A[a++];return O}
for(;++i<4;)L=(L<<8|A[a++])>>>0;
//read word lengths
let D=[];
for(let G=Uint16Array.from(Array(60),a=>2048),a=i=s=0;i<kn;){
let l=0,b=0;
for(;!db(a,G,5);)a++;
for(l=1<<a;a;)l|=db(++b+29,G,6)<<--a;
s+=D[i++]=l
}
//build dict
let cs=i=0,E,F=new Uint32Array(c=kn>256?kn+1:257),C=new Uint32Array(c+1);
for(gets(F,256,s);i<s;)C[cs+1]=i+=F[cs++];
for(E=F.slice(i=0);s;){
for(let l=D[i],p=0,W=D[i++]=new Uint8Array(l);p<l;){
for(;R<16777216;R*=256)L=(L<<8|A[a++])>>>0;
R=R/s>>>0;
let b=0,c=cs-1,d,e=L/R>>>0;
for(;b<c;C[d+1]>e?c=d:b=++d)d=b+c>>1;
L-=R*C[W[p++]=b],R*=F[b];
if(!--E[b])for(b=s=0;b<cs;)C[b+1]=s+=F[b]=E[b++]
}
}
//extract model
for(let G=Uint16Array.from(Array(60),a=>2048),a=0,p=kn;p;){
let l=db(29,G,4),b=0;
if(l){
for(;!db(a,G,5);)a++;
for(l=1<<a;a;)l|=db(++b+29,G,6)<<--a
}F[--p]=l
}
for(i=kn,s=0;i;)s=F[i]=s+F[--i];
for(i=s=0;i<kn;)C[i+1]=s+=E[i]=F[i++]=F[i];
//extract all
for(i=0;s;){
for(;R<16777216;R*=256)L=(L<<8|A[a++])>>>0;
R=R/s>>>0;
let b=0,c=kn-1,d,e=L/R>>>0;
for(;b<c;C[d+1]>e?c=d:b=++d)d=b+c>>1;
L-=R*C[b],R*=F[b];
if(!--E[b])for(c=s=0;c<kn;)C[c+1]=s+=F[c]=E[c++];
O.push.apply(O,D[b])
}return O
}
ちなみにentropyTableは気持悪い形式で圧縮してあります。文字符号はShiftJISとしておきます。
それではcodepen召喚
See the Pen Data compression by Word eXtraction method by xezz (@xezz) on CodePen.
本家はとっくの昔に消滅したのでここに晒しておきます