testlag.m

HMMBOX, version 3.2, William Penny, Imperial College, Feb 1999 Matlab toolbox for Hidden Markov Mod

P=rand(2,4); 
for i=1:4, P(:,i)=P(:,i)./sum(P(:,i)); end;
for i=1:4, P2(:,i)=reshape(P(:,i)*P(:,i)',4,1); end;
hmm1.P=P2;
hmm1.K=4;
N=4096;
hmm1.Pi=ones(1,hmm1.K)./hmm1.K;
hmm1.obsmodel='Gauss';
hmm1.state(1).Mu=-6;hmm1.state(2).Mu=-3;hmm1.state(3).Mu=3;hmm1.state(4).Mu=6;
hmm1.state(1).Mu=-6;hmm1.state(2).Mu=-3;hmm1.state(3).Mu=3;hmm1.state(4).Mu=6;
hmm1.state(1).Cov=1.5;hmm1.state(2).Cov=1.5;hmm1.state(3).Cov=1.5;hmm1.state(4).Cov=1.5;
simdata1=hmmsim(hmm1,N);
hmm2.Pi=[1/2 1/2];
hmm2.P=[3/4 1/2; 1/4 1/2];
hmm2.K=2;
hmm2.obsmodel='Gauss';
hmm2.state(1).Mu=-5;hmm2.state(2).Mu=0;
hmm2.state(1).Cov=2;hmm2.state(2).Cov=1;
simdata2=hmmsim(hmm2,N);
cd ../complex
X1=simdata1.Xclass;
X2=simdata2.Xclass;
U1=embed(X1,20,1);
S1=svd(U1,0);
U2=embed(X2,20,1);
S2=svd(U2,0);
plot(1:20,S1,1:20,S2)
U1=embed(X1,20,2);
S1=svd(U1,0);
plot(1:20,S1,1:20,S2)
U1=embed(X1,20,3);
S1=svd(U1,0);
plot(1:20,S1,1:20,S2)
U1=embed(X1,20,4);
S1=svd(U1,0);
plot(1:20,S1,1:20,S2)


hmm1.state(1).Mu=-5;hmm1.state(2).Mu=0;hmm1.state(3).Mu=-5; ...
    hmm1.state(4).Mu=0;
cd ../hmmbox_3_2

simdata1=hmmsim(hmm1,N);
hmmt1.K=2;
hmmt1=hmminit(simdata1.Xseries,hmmt1,'diag');
hmmt1.train.cyc=30;
hmmt1.obsmodel='Gauss';
hmmt1.train.obsupdate=ones(1,hmmt1.K);  
hmmt1.train.init=1;       
hmmt1=hmmtrain(simdata1.Xseries,N,hmmt1);
[block1]=hmmdecode(simdata1.Xseries,N,hmmt1);

hmmt2.K=2;
hmmt2=hmminit(simdata2.Xseries,hmmt2,'diag');
hmmt2.train.cyc=30;
hmmt2.obsmodel='Gauss';
hmmt2.train.obsupdate=ones(1,hmmt2.K);  
hmmt2.train.init=1;       
hmmt2=hmmtrain(simdata2.Xseries,N,hmmt2);
[block2]=hmmdecode(simdata2.Xseries,N,hmmt2);

cd ../complex

X1=block1.q_star;
X2=block2.q_star;
U1=embed(X1,20,1);
S1=svd(U1,0);
U2=embed(X2,20,1);
S2=svd(U2,0);
plot(1:20,S1,1:20,S2)
U1=embed(X1,20,2);
S1=svd(U1,0);
plot(1:20,S1,1:20,S2)
U1=embed(X1,20,3);
S1=svd(U1,0);
plot(1:20,S1,1:20,S2)
U1=embed(X1,20,4);
S1=svd(U1,0);
plot(1:20,S1,1:20,S2)
U2=embed(X2,20,4);
S2=svd(U2,0);
plot(1:20,S1,1:20,S2)

cd ../gmm
[gmm1,lik] = gmmem (struct('timeseries',simdata1.Xseries),2);
[gmm2,lik] = gmmem (struct('timeseries',simdata2.Xseries),2);


cd ../complex

X1=gmm1.pjgx;
X2=gmm2.pjgx;
U1=embed(X1,20,1);
S1=svd(U1,0);
U2=embed(X2,20,1);
S2=svd(U2,0);
plot(1:20,S1,1:20,S2)
U1=embed(X1,20,2);
S1=svd(U1,0);
plot(1:20,S1,1:20,S2)
U1=embed(X1,20,3);
S1=svd(U1,0);
plot(1:20,S1,1:20,S2)
U1=embed(X1,20,4);
S1=svd(U1,0);
plot(1:20,S1,1:20,S2)
U2=embed(X2,20,4);
S2=svd(U2,0);
plot(1:20,S1,1:20,S2)
U1=embed(X1,10,1);
S1=svd(U1,0);
U2=embed(X2,10,1);
S2=svd(U2,0);
plot(1:10,S1,1:10,S2)