%
%tested 1 in March 26.2002 P/181, Graphical models in applied multivariate statistics
%test mathematic marks, 1:mech, 2:vect, 3:alg, 4:anal, 5:stat
clear all
close all
testdata=[77 82 67 67 81;...
63 78 80 70 81;...
75 73 71 66 81;...
55 72 63 70 68;...
63 63 65 70 63;...
53 61 72 64 73;...
51 67 65 65 68;...
59 70 68 62 56;...
62 60 58 62 70;...
64 72 60 62 45;...
52 64 60 63 54;...
55 67 59 62 44;...
50 50 64 55 63;...
65 63 58 56 37;...
31 55 60 57 73;...
60 64 56 54 40;...
44 69 53 53 53;...
42 69 61 55 45;...
62 46 61 57 45;...
31 49 62 63 62;...
44 61 52 62 46;...
49 41 61 49 64;...
12 58 61 63 67;...
49 53 49 62 47;...
54 49 56 47 53;...
54 53 46 59 44;...
44 56 55 61 36;...
18 44 50 57 81;...
46 52 65 50 35;...
32 45 49 57 64;...
30 69 50 52 45;...
46 49 53 59 37;...
40 27 54 61 61;...
31 42 48 54 68;...
36 59 51 45 51;...
56 40 56 54 35;...
46 56 57 49 32;...
45 42 55 56 40;...
42 60 54 49 33;...
40 63 53 54 25;...
23 55 59 53 44;...
48 48 49 51 37;...
41 63 49 46 34;...
46 52 53 41 40;...
46 61 46 38 41;...
40 57 51 52 31;...
49 49 45 48 39;...
22 58 53 56 41;...
35 60 47 54 33;...
48 56 49 42 32;...
31 57 50 54 34;...
17 53 57 43 51;...
49 57 47 39 26;...
59 50 47 15 46;...
37 56 49 28 45;...
40 43 48 21 61;...
35 35 41 51 50;...
38 44 54 47 24;...
43 43 38 34 49;...
39 46 46 32 43;...
62 44 36 22 42;...
48 38 41 44 33;...
34 42 50 47 29;...
18 51 40 56 30;...
35 36 46 48 29;...
59 53 37 22 19;...
41 41 43 30 33;...
31 52 37 27 40;...
17 51 52 35 31;...
34 30 50 47 36;...
46 40 47 29 17;...
10 46 36 47 39;...
46 37 45 15 30;...
30 34 43 46 18;...
13 51 50 25 31;...
49 50 38 23 9;...
18 32 31 45 40;...
8 42 48 26 40;...
23 38 36 48 15;...
30 24 43 33 25;...
3 9 51 47 40;...
7 51 43 17 22;...
15 40 43 23 18;...
15 38 39 28 17;...
5 30 44 36 18;...
12 30 32 35 21;...
5 26 15 20 20;...
0 40 21 9 14 ];
testdata=som_normalize(testdata,'var');
%old_var=[ 302.29 125.78 100.43 105.07 116.07; ...
% 125.78 170.88 84.19 93.6 97.89;...
% 100.43 84.19 111.6 110.84 120.49; ...
% 105.07 93.6 110.84 217.88 153.77; ...
% 116.07 97.89 120.49 153.77 294.37]
old_var=cov(testdata)
N=88 % number of observations
q=5 % number of variables
cutoff=0.05;
[old_parcorrf,delet_edge,record_edge]=MGraph_GaussModel(old_var,N,q,cutoff)
spareM=(abs(old_parcorrf)>0.0000001 & abs(old_parcorrf)~=1);
spareM=sparse(spareM);
draw_graph(spareM);
%testd 2 Introduction to Graphical Modelling P/44
%test DigoxinClearance
%it is correct when device set to d=abs(N*log(det(old_var)/det(fit_var)))
%and normalize the raw data
close all
clear all
data=[19.5 17.5 0.74; 24.7 34.8 0.43;26.5 11.4 0.11; 31.1 29.3 1.48;...
31.3 13.9 0.97; 31.8 31.6 1.12; 34.1 20.7 1.77; 36.6 34.1 0.7;...
42.4 25 0.93; 42.8 47.4 2.5; 44.2 31.8 0.89; 49.7 36.1 0.52; ...
51.3 22.7 0.33; 55 30.7 0.8; 55.9 42.5 1.02; 61.2 42.4 0.56; 63.1 61.1 0.93; ...
63.7 38.2 0.44; 66.8 37.5 0.5; 72.4 50.1 0.97; 80.9 50.2 1.02; 82 50 0.95;...
82.7 31.8 0.76; 87.9 55.4 1.06; 101.5 110.6 1.38; 105 114.4 1.85; ...
110.5 69.3 2.25 ; 114.2 84.8 1.76; 117.8 63.9 1.6; 122.6 76.1 0.88; ...
127.9 112.8 1.7; 135.6 82.2 0.98 ; 136 46.8 0.94; 153.5 137.7 1.76; 201.1 76.1 0.87]
data=som_normalize(data,'var');
N=35 %number of observations
q=3 %number of dimensional
cutoff=0.05
[old_parcorrf,delet_edge,record_edge]=MGraph_GaussModel(cov(data),N,q,cutoff)
%make a sparse matrix to plot the network
spareM=(abs(old_parcorrf)>0.0000001 & abs(old_parcorrf)~=1);
spareM=sparse(spareM);
draw_graph(spareM);
%needmore test 3
%the reason is ,it should be decomposal model
%need one more function to detect decomposable model !!
%
%anxity and anger
%test not match the results of MIM, both device choic are not match with MIM
%normalize of raw data also not match MIM
%maybe something wrong with the spares matrise!!!!
%need check with introduction graphical model...
%cont WXYZ
%label W "Anxiety st" X "Anger st" Y "Anxiety tr" Z "Anger tr"
%sread WXYZ
%bug at here, must adjust the degree of freedom for sparse(zero) entry
clear all
close all
%count=[684
%mean=[18.8744, 15.2265, 21.2019, 23.4217 ;...
cov_data=[37.1926 24.9311 21.6056 15.6907 ; ...
24.9311 44.8472 17.8072 21.8565 ;...
21.6056 17.8072 32.2462 18.3523;...
15.6907 21.8565 18.3523 43.1191 ];
%data=som_normalize(data,'var');
N=684 %number of observations
q=4 %number of dimensional
cutoff=0.05
[old_parcorrf,delet_edge,record_edge]=MGraph_GaussModel(cov_data,N,q,cutoff)
spareM=(abs(old_parcorrf)>0.0000001 & abs(old_parcorrf)~=1);
spareM=sparse(spareM);
draw_graph(spareM);
%needmore test 4
%test failed with raw data but it works for the normalized data,
%maybe need checked the device and df....
%that means if the data in same scale dont need normalize but differ scales need
%normalize! Need check with introduction graphical model
%Bone mineartest
%cont UVWXYZ
%label U "Men age" V "BMI" W "Alk" X "BMC-arm" Y "BMD-sp" Z "BMC-sp"
%read UVWXYZ
clear all
close all
data=[ 59.000 0.002 4.963 1.110 1.068 4.740;...
86.000 0.003 5.333 1.370 1.039 4.480;...
116.000 0.002 5.247 1.320 0 0;...
129.000 0.002 5.323 1.500 0 0;...
59.000 0.003 5.283 1.090 0.930 4.390;...
49.000 0.002 4.977 0.980 0.785 3.420;...
75.000 0.003 4.963 1.160 1.217 5.140;...
51.000 0.002 5.017 1.120 0.762 3.250;...
62.000 0.002 4.984 1.210 0.892 3.970;...
123.000 0.002 4.868 1.560 1.007 5.260;...
91.000 0.002 5.124 1.500 0 0;...
53.000 0.002 4.990 1.130 0.907 3.770;...
70.000 0.002 4.927 1.310 0.711 3.250;...
104.000 0.002 4.852 1.000 0.934 3.840;...
72.000 0.002 4.700 1.490 1.248 5.380;...
84.000 0.003 5.182 1.100 1.037 4.580;...
101.000 0.003 4.990 1.100 0.881 4.420;...
85.000 0.002 4.920 1.540 0.968 4.340;...
80.000 0.003 5.220 0.930 1.069 4.690;...
51.000 0.002 4.868 1.050 0.875 3.800;...
85.000 0.002 5.004 1.290 0.752 3.310;...
71.000 0.003 5.204 0.870 0.834 3.680;...
105.000 0.002 4.883 1.130 1.109 5.170;...
101.000 0.003 4.949 1.130 0.972 4.120;...
94.000 0.003 5.118 1.180 0.871 4.200;...
76.000 0.002 4.691 1.150 1.069 5.330;...
64.000 0.003 5.193 1.260 0.938 4.420;...
113.000 0.003 5.293 1.050 0.991 4.710;...
124.000 0.002 5.357 0.960 0.885 3.570;...
103.000 0.002 4.828 1.060 0.949 4.090;...
106.000 0.002 4.745 1.270 0.886 3.710;...
116.000 0.002 5.293 0.940 0.936 4.200;...
87.000 0.002 5.176 0.960 0.810 4.190;...
133.000 0.003 4.934 1.260 1.135 5.200;...
90.000 0.002 4.762 0.940 0.974 4.310;...
75.000 0.003 4.812 1.070 0.987 4.370;...
123.000 0.003 5.263 0.890 0.716 3.430;...
82.000 0.002 5.568 1.220 0.693 3.150;...
134.000 0.002 5.236 1.010 0.875 3.780;...
112.000 0.002 5.347 1.020 1.019 5.010;...
123.000 0.003 4.771 1.380 1.041 4.870;...
134.000 0.002 5.124 1.120 0.693 2.990;...
130.000 0.002 5.226 0.880 0.948 4.230;...
114.000 0.002 5.153 1.200 0.925 4.060;...
97.000 0.002 4.970 1.230 0.805 3.520;...
54.000 0.002 5.100 1.160 0.818 3.650;...
58.000 0.002 4.852 1.590 1.151 5.060;...
82.000 0.003 5.081 1.500 1.133 4.960;...
140.000 0.002 5.198 1.170 0.763 3.600;...
75.000 0.002 4.691 1.470 0 0;...
79.000 0.003 4.727 1.270 1.036 5.080;...
111.000 0.002 5.142 1.040 0.771 3.420;...
108.000 0.003 5.587 1.120 0.906 4.140;...
136.000 0.002 5.659 1.130 0.685 3.140;...
51.000 0.002 4.956 1.150 0 0;...
125.000 0.003 4.997 1.200 0 0;...
51.000 0.002 5.170 1.430 0.886 4.070;...
88.000 0.002 5.204 0.940 0.790 3.410;...
87.000 0.002 5.124 1.110 1.093 5.090;...
66.000 0.002 4.913 1.060 0.915 3.780;...
59.000 0.002 5.278 1.300 0.979 4.530;...
73.000 0.002 5.159 1.180 0.867 4.240;...
105.000 0.003 5.328 1.050 0.853 3.660;...
61.000 0.003 4.852 1.070 0.901 3.870;...
99.000 0.002 5.011 1.280 0.916 4.070;...
53.000 0.002 5.293 1.400 0.999 4.690;...
141.000 0.003 5.165 1.090 0.908 4.270;...
69.000 0.002 5.533 1.170 0.784 3.550;...
103.000 0.002 5.501 0.960 0.624 2.620;...
49.000 0.002 4.963 0.880 0.867 4.100;...
123.000 0.002 5.124 1.100 0.927 4.300;...
83.000 0.003 5.204 1.150 0.972 4.280;...
84.000 0.002 5.342 1.290 0.885 4.200;...
104.000 0.002 4.949 1.210 0.929 3.960;...
126.000 0.002 4.787 1.330 0.865 4.090;...
87.000 0.002 4.990 1.140 0.895 4.470;...
126.000 0.002 5.298 1.520 1.119 5.450;...
102.000 0.002 5.043 1.820 1.326 6.250;...
79.000 0.002 4.804 1.050 0.954 4.050;...
52.000 0.002 4.913 1.240 1.005 4.610;...
135.000 0.002 5.142 1.220 0.831 3.910;...
71.000 0.002 4.898 0.900 0.659 2.870;...
112.000 0.002 5.257 1.520 0.817 3.710;...
78.000 0.002 4.890 0.990 0.801 3.660;...
123.000 0.002 5.069 1.030 0.807 3.310;...
81.000 0.002 5.303 0.890 0.933 4.230;...
77.000 0.003 5.182 1.410 1.014 4.530;...
134.000 0.002 4.997 1.030 0.964 4.370;...
128.000 0.002 5.165 0.900 0.837 3.950;...
116.000 0.002 5.313 1.400 0.995 4.400;...
108.000 0.002 4.828 1.290 0.915 3.870;...
95.000 0.002 5.069 1.480 1.147 5.350;...
70.000 0.002 4.913 1.040 0.823 3.690;...
80.000 0.003 5.247 0.990 0.975 4.560;...
121.000 0.002 5.011 1.290 0.927 3.760;...
43.000 0.003 4.875 1.030 0.927 4.320;...
77.000 0.003 5.142 1.020 0.841 4.070;...
127.000 0.003 5.182 1.310 1.151 5.730;...
85.000 0.003 5.352 0.980 0.927 4.130;...
105.000 0.003 5.352 0.990 0.994 4.180;...
119.000 0.002 4.920 1.320 1.194 5.090;...
79.000 0.002 5.030 1.300 1.271 5.720;...
68.000 0.002 5.361 1.120 0 0;...
134.000 0.003 4.984 0.930 0 0;...
93.000 0.002 4.883 1.350 1.005 4.010;...
82.000 0.002 5.204 1.260 0.946 4.390;...
136.000 0.003 5.398 1.020 0.756 3.370;...
78.000 0.002 5.094 1.050 0.803 3.670;...
55.000 0.003 5.220 1.310 1.173 5.130;...
104.000 0.002 5.050 1.170 0.868 3.950;...
49.000 0.002 5.403 0 0.857 4.020;...
132.000 0.002 5.476 0.980 0.932 4.260;...
75.000 0.002 4.927 0.970 0.787 3.550;...
136.000 0.003 5.112 0.940 0.842 4.080;...
86.000 0.003 5.323 1.020 0 0;...
133.000 0.002 5.024 1.190 1.064 4.820;...
106.000 0.002 5.263 1.010 0.821 3.640;...
65.000 0.002 5.323 1.360 1.114 5.510;...
62.000 0.002 5.165 1.240 0.889 3.730;...
124.000 0.003 5.565 0.820 0.695 3.190;...
128.000 0.002 4.942 1.060 1.024 4.750;...
130.000 0.002 5.130 1.320 0 0;...
64.000 0.002 5.198 1.450 0.990 4.510;...
64.000 0.002 4.913 1.460 1.252 5.830;...
102.000 0.003 5.063 1.050 1.029 4.570;...
64.000 0.002 4.736 1.120 0.900 4.440;...
78.000 0.002 4.883 1.300 0.995 4.590;...
67.000 0.002 5.170 1.400 0.960 4.600;...
99.000 0.002 4.927 1.040 0.906 3.910;...
53.000 0.002 4.820 1.070 1.054 4.870;...
79.000 0.002 4.905 1.020 0.920 3.970;...
59.000 0.003 5.252 1.460 1.081 5.870;...
93.000 0.002 5.130 1.150 0.919 3.930;...
53.000 0.002 5.371 1.000 0.781 3.780;...
87.000 0.002 5.318 1.330 0.971 4.420;...
130.000 0.002 4.585 1.200 0.859 4.260;...
111.000 0.003 5.142 0.980 0.878 3.780;...
130.000 0.002 5.273 1.310 1.031 5.060;...
139.000 0.003 5.242 0.990 0.838 3.800;...
106.000 0.002 5.398 0.840 0.832 3.730;...
131.000 0.003 5.635 1.130 0.799 3.740;...
81.000 0.002 5.043 1.360 0.962 4.310;...
130.000 0.002 5.288 1.210 0.947 4.310;...
71.000 0.002 4.828 1.340 1.006 5.360;...
135.000 0.003 5.050 1.130 0.703 3.340;...
64.000 0.002 4.852 1.430 0.959 4.110;...
63.000 0.003 5.293 1.210 0.996 4.480;...
78.000 0.002 5.236 1.110 0.815 3.650;...
59.000 0.002 4.942 1.020 0.805 3.800;...
103.000 0.003 5.142 1.000 0.788 3.380 ];
data=som_normalize(data,'var') ; %normalize data
N=150 %number of observations
q=6 %number of dimensional
cutoff=0.05
[old_parcorrf,delet_edge,record_edge]=MGraph_GaussModel(cov(data),N,q,cutoff)
spareM=(abs(old_parcorrf)>0.0000001 & abs(old_parcorrf)~=1);
spareM=sparse(spareM);
draw_graph(spareM);
%test 5
%test Diagon clearance
%when device choose d=abs(N*log(det(old_var)/det(fit_var))) and normalize the raw data
%it match the MIM results
% Data from Altman, p. 323
% x=Creatinine clearance, y=Digoxin clearance, z=urine flow
%cont xyz
%read xyz
clear all
close all
data=[ 19.5 17.5 0.74 ;...
24.7 34.8 0.43;...
26.5 11.4 0.11;...
31.1 29.3 1.48;...
31.3 13.9 0.97;...
31.8 31.6 1.12;...
34.1 20.7 1.77;...
36.6 34.1 0.70;...
42.4 25.0 0.93;...
42.8 47.4 2.50;...
44.2 31.8 0.89;...
49.7 36.1 0.52;...
51.3 22.7 0.33;...
55.0 30.7 0.80;...
55.9 42.5 1.02;...
61.2 42.4 0.56;...
63.1 61.1 0.93;...
63.7 38.2 0.44;...
66.8 37.5 0.50;...
72.4 50.1 0.97;...
80.9 50.2 1.02;...
82.0 50.0 0.95;...
82.7 31.8 0.76;...
87.9 55.4 1.06;...
101.5 110.6 1.38;...
105.0 114.4 1.85;...
110.5 69.3 2.25;...
114.2 84.8 1.76;...
117.8 63.9 1.60;...
122.6 76.1 0.88;...
127.9 112.8 1.70;...
135.6 82.2 0.98;...
136.0 46.8 0.94;...
153.5 137.7 1.76;...
201.1 76.1 0.87];
data=som_normalize(data,'var') ; %normalize data
N=35 %number of observations
q=3 %number of dimensional
cutoff=0.05
[old_parcorrf,delet_edge,record_edge]=MGraph_GaussModel(cov(data),N,q,cutoff)
spareM=(abs(old_parcorrf)>0.0000001 & abs(old_parcorrf)~=1);
spareM=sparse(spareM);
draw_graph(spareM);
%test 6
%needmore test fret's head
%need dectec the decomposable model!!
% Fret's heads data.
%P/255 Graphical models in applied multivariant statistical analysis
%both device choose of test did not match with MIM!
%the device value was not correct !! nov 12.2002
%cont wxyz
%label w "1st Length " x "1st Breadth" y "2nd Length" z "2nd Breadth"
%read wxyz
clear all
close all
data=[191 155 179 145;...
195 149 201 152;...
181 148 185 149;...
183 153 188 149;...
176 144 171 142;...
208 157 192 152;...
189 150 190 149;...
197 159 189 152;...
188 152 197 159;...
192 150 187 151;...
179 158 186 148;...
183 147 174 147;...
174 150 185 152;...
190 159 195 157;...
188 151 187 158;...
163 137 161 130;...
195 155 183 158;...
186 153 173 148;...
181 145 182 146;...
175 140 165 137;...
192 154 185 152;...
174 143 178 147;...
176 139 176 143;...
197 167 200 158;...
190 163 187 150 ]
%data=som_normalize(data,'var');
N=25 %number of observations
q=4 %number of dimensional
cutoff=0.05
[old_parcorrf,delet_edge,record_edge]=MGraph_GaussModel(cov(data),N,q,cutoff)
spareM=(abs(old_parcorrf)>0.0000001 & abs(old_parcorrf)~=1);
spareM=sparse(spareM);
draw_graph(spareM);
%needmore test 8 IRIS data
%results match MIM but the device value are different?? Nov12.2002
% Data that should need no introduction.
% Three species of iris: Iris Setosa (a=1), Iris versicolour (a=2),
% Iris virginia (a=3).
%cont wxyz
%label w "Sepal length" x "Sepal width" y "Petal length" z "Petal width"
%read wxyz
close all
clear all
data=[5.1 3.5 1.4 0.2;...
4.9 3.0 1.4 0.2;...
4.7 3.2 1.3 0.2;...
4.6 3.1 1.5 0.2;...
5.0 3.6 1.4 0.2;...
5.4 3.9 1.7 0.4;...
4.6 3.4 1.4 0.3;...
5.0 3.4 1.5 0.2;...
4.4 2.9 1.4 0.2;...
4.9 3.1 1.5 0.1;...
5.4 3.7 1.5 0.2;...
4.8 3.4 1.6 0.2;...
4.8 3.0 1.4 0.1;...
4.3 3.0 1.1 0.1;...
4.8 4.0 1.2 0.2;...
5.7 4.4 1.5 0.4;...
5.4 3.9 1.3 0.4;...
5.1 3.5 1.4 0.3;...
5.7 3.8 1.7 0.3;...
5.1 3.8 1.5 0.3;...
5.4 3.4 1.7 0.2;...
5.1 3.7 1.5 0.4;...
4.6 3.6 1.0 0.2;...
5.1 3.3 1.7 0.5;...
4.8 3.4 1.9 0.2;...
5.0 3.0 1.6 0.2;...
5.0 3.4 1.6 0.4;...
5.2 3.5 1.5 0.2;...
5.2 3.4 1.4 0.2;...
4.7 3.2 1.6 0.2;...
4.8 3.1 1.6 0.2;...
5.4 3.4 1.5 0.4;...
5.2 4.1 1.5 0.1;...
5.5 4.2 1.4 0.2;...
4.9 3.1 1.5 0.2;...
5.0 3.2 1.2 0.2;...
5.5 3.5 1.3 0.2;...
4.9 3.6 1.4 0.1;...
4.4 3.0 1.3 0.2;...
5.1 3.4 1.5 0.2;...
5.0 3.5 1.3 0.3;...
4.5 2.3 1.3 0.3;...
4.4 3.2 1.3 0.2;...
5.0 3.5 1.6 0.6;...
5.1 3.8 1.9 0.4;...
4.8 3.0 1.4 0.3;...
5.1 3.8 1.6 0.2;...
4.6 3.2 1.4 0.2;...
5.3 3.7 1.5 0.2;...
5.0 3.3 1.4 0.2;...
7.0 3.2 4.7 1.4;...
6.4 3.2 4.5 1.5;...
6.9 3.1 4.9 1.5;...
5.5 2.3 4.0 1.3;...
6.5 2.8 4.6 1.5;...
5.7 2.8 4.5 1.3;...
6.3 3.3 4.7 1.6;...
4.9 2.4 3.3 1.0;...
6.6 2.9 4.6 1.3;...
5.2 2.7 3.9 1.4;...
5.0 2.0 3.5 1.0;...
5.9 3.0 4.2 1.5;...
6.0 2.2 4.0 1.0;...
6.1 2.9 4.7 1.4;...
5.6 2.9 3.6 1.3;...
6.7 3.1 4.4 1.4;...
5.6 3.0 4.5 1.5;...
5.8 2.7 4.1 1.0;...
6.2 2.2 4.5 1.5;...
5.6 2.5 3.9 1.1;...
5.9 3.2 4.8 1.8;...
6.1 2.8 4.0 1.3;...
6.3 2.5 4.9 1.5;...
6.1 2.8 4.7 1.2;...
6.4 2.9 4.3 1.3;...
6.6 3.0 4.4 1.4;...
6.8 2.8 4.8 1.4;...
6.7 3.0 5.0 1.7;...
6.0 2.9 4.5 1.5;...
5.7 2.6 3.5 1.0;...
5.5 2.4 3.8 1.1;...
5.5 2.4 3.7 1.0;...
5.8 2.7 3.9 1.2;...
6.0 2.7 5.1 1.6;...
5.4 3.0 4.5 1.5;...
6.0 3.4 4.5 1.6;...
6.7 3.1 4.7 1.5;...
6.3 2.3 4.4 1.3;...
5.6 3.0 4.1 1.3;...
5.5 2.5 4.0 1.3;...
5.5 2.6 4.4 1.2;...
6.1 3.0 4.6 1.4;...
5.8 2.6 4.0 1.2;...
5.0 2.3 3.3 1.0;...
5.6 2.7 4.2 1.3;...
5.7 3.0 4.2 1.2;...
5.7 2.9 4.2 1.3;...
6.2 2.9 4.3 1.3;...
5.1 2.5 3.0 1.1;...
5.7 2.8 4.1 1.3;...
6.3 3.3 6.0 2.5;...
5.8 2.7 5.1 1.9;...
7.1 3.0 5.9 2.1;...
6.3 2.9 5.6 1.8;...
6.5 3.0 5.8 2.2;...
7.6 3.0 6.6 2.1;...
4.9 2.5 4.5 1.7;...
7.3 2.9 6.3 1.8;...
6.7 2.5 5.8 1.8;...
7.2 3.6 6.1 2.5;...
6.5 3.2 5.1 2.0;...
6.4 2.7 5.3 1.9;...
6.8 3.0 5.5 2.1;...
5.7 2.5 5.0 2.0;...
5.8 2.8 5.1 2.4;...
6.4 3.2 5.3 2.3;...
6.5 3.0 5.5 1.8;...
7.7 3.8 6.7 2.2;...
7.7 2.6 6.9 2.3;...
6.0 2.2 5.0 1.5;...
6.9 3.2 5.7 2.3;...
5.6 2.8 4.9 2.0;...
7.7 2.8 6.7 2.0;...
6.3 2.7 4.9 1.8;...
6.7 3.3 5.7 2.1;...
7.2 3.2 6.0 1.8;...
6.2 2.8 4.8 1.8;...
6.1 3.0 4.9 1.8;...
6.4 2.8 5.6 2.1;...
7.2 3.0 5.8 1.6;...
7.4 2.8 6.1 1.9;...
7.9 3.8 6.4 2.0;...
6.4 2.8 5.6 2.2;...
6.3 2.8 5.1 1.5;...
6.1 2.6 5.6 1.4;...
7.7 3.0 6.1 2.3;...
6.3 3.4 5.6 2.4;...
6.4 3.1 5.5 1.8;...
6.0 3.0 4.8 1.8;...
6.9 3.1 5.4 2.1;...
6.7 3.1 5.6 2.4;...
6.9 3.1 5.1 2.3;...
5.8 2.7 5.1 1.9;...
6.8 3.2 5.9 2.3;...
6.7 3.3 5.7 2.5;...
6.7 3.0 5.2 2.3;...
6.3 2.5 5.0 1.9;...
6.5 3.0 5.2 2.0;...
6.2 3.4 5.4 2.3;...
5.9 3.0 5.1 1.8];
data=som_normalize(data,'var');
N=150 %number of observations
q=4 %number of dimensional
cutoff=0.05
[old_parcorrf,delet_edge,record_edge]=MGraph_GaussModel(cov(data),N,q,cutoff)
spareM=(abs(old_parcorrf)>0.0000001 & abs(old_parcorrf)~=1);
spareM=sparse(spareM);
draw_graph(spareM);
%test 9 mice data
%when device choose d=abs(N*log(det(old_parcorrf)/det(fit_parcorrf)))
%it is match the MIM resuls
%But after apply normalize the raw data to mean =0, var=1 it match with MIM well
% Mice Data.
% Source: Morrison (1976). Multivariate Statistical Methods, McGraw-Hill.
%
%cont xyz
%read xyz
%test matchs MIM when change device to d=abs(N*log(det(old_parcorrf)/det(fit_parcorrf)))
clear all
close all
data=[1.21 0.61 0.70;...
0.92 0.43 0.71;...
0.80 0.43 0.71;...
0.85 0.48 0.68;...
0.98 0.42 0.71;...
1.15 0.52 0.72;...
1.10 0.50 0.75;...
1.02 0.53 0.70;...
1.18 0.45 0.70;...
1.09 0.40 0.69;...
1.40 0.50 0.71;...
1.17 0.39 0.69;...
1.23 0.44 0.70;...
1.19 0.37 0.72;...
1.38 0.42 0.71;...
1.17 0.45 0.70;...
1.31 0.41 0.70;...
1.30 0.47 0.67;...
1.22 0.29 0.68;...
1.00 0.30 0.70;...
1.12 0.27 0.72;...
1.09 0.35 0.73]
data=som_normalize(data,'var');
N=22 %number of observations
q=3 %number of dimensional
cutoff=0.05
[old_parcorrf,delet_edge,record_edge]=MGraph_GaussModel(cov(data),N,q,cutoff)
spareM=(abs(old_parcorrf)>0.0000001 & abs(old_parcorrf)~=1);
spareM=sparse(spareM);
draw_graph(spareM);
%test latian square
clear all
load latinsquare.dat
red_1T1=latinsquare(:,1);
red_1C1=latinsquare(:,2);
green_2T2=latinsquare(:,3);
green_2c2=latinsquare(:,4);
data=log([latinsquare(:,1:4)]);
data=som_normalize(data,'var');
N=336;
q=4;
cutoff=0.05
[old_parcorrf,delet_edge,record_edge]=MGraph_GaussModel(cov(data),N,q,cutoff)
spareM=(abs(old_parcorrf)>0.0000001 & abs(old_parcorrf)~=1);
spareM=sparse(spareM);
draw_graph(spareM);
%test reference design
%
%test latian square
clear all
load refdesign.dat
data=([refdesign(:,1:4)]);
data=som_normalize(data,'var');
N=1905
q=4;
cutoff=0.05
[old_parcorrf,delet_edge,record_edge]=MGraph_GaussModel(cov(data),N,q,cutoff)
spareM=(abs(old_parcorrf)>0.0000001 & abs(old_parcorrf)~=1);
spareM=sparse(spareM);
draw_graph(spareM);
%test niesh
%failed
clear all
load niehs.dat
data=(niehs);
data=som_normalize(data,'var');
N=1907;
q=12;
cutoff=0.05
[old_parcorrf,delet_edge,record_edge]=MGraph_GaussModel(cov(data),N,q,cutoff)
spareM=(abs(old_parcorrf)>0.0000001 & abs(old_parcorrf)~=1);
spareM=sparse(spareM);
[x,y,h]=draw_graph(spareM);
%test dnr_jmaerry data
%mouse size with expression level
%
clear all
close all
filename='Merger23_filesID.txt';
[tdata, varnames, casenames] = tblread(filename,'\t') ;
data=log2(tdata);
%1343 genes, 23 samples
t_data=data(:,[1,2,3,4,5]); %test mouse size with expression level
t_data=som_normalize(t_data,'var');
N=1343;
q=5;
cutoff=0.01
[old_parcorrf,delet_edge,record_edge]=MGraph_GaussModel(cov(t_data),N,q,cutoff)
spareM=(abs(old_parcorrf)>0.0000001 & abs(old_parcorrf)~=1);
spareM=sparse(spareM);
[x,y,h]=draw_graph(spareM);
%test PKC path way from science 2000
clear all
close all
filename='Pkc_pathData.dat';
[tdata, varnames, casenames] = tblread(filename,'\t') ;
%t_data=tdata([1,4,5:8,10,11,13],:)';
t_data=tdata';
t_data=som_normalize(t_data,'var');
N=46;
q=13;
cutoff=0.05
[old_parcorrf,delet_edge,record_edge]=MGraph_GaussModel(cov(t_data),N,q,cutoff)
spareM=(abs(old_parcorrf)>0.0000001 & abs(old_parcorrf)~=1);
spareM=sparse(spareM);
genenames=strrep(cellstr( deblank(casenames) ),'"','');
genenames(1)={'wsc1'};
genenames(4)={'mid2'};
genenames(6)={'pkc1'};
genenames(7)={'bck1'};
genenames(8)={'mkk1'};
genenames(9)={'mkk2'};
genenames(10)={'mpk1'};
genenames(11)={'swi4'};
genenames(12)={'swi6'};
[x,y,h]=draw_graph(spareM,lower(genenames));
%test Pheromone path way from science 2000
clear all
close all
filename='pheromone_pathway.dat';
[tdata, varnames, casenames] = tblread(filename,'\t') ;
%t_data=tdata([1,4:11,13:14],:);
t_data=tdata;
t_data=som_normalize(t_data','var');
N=46;
q=15;
cutoff=0.05
[old_parcorrf,delet_edge,record_edge]=MGraph_GaussModel(cov(t_data),N,q,cutoff)
spareM=(abs(old_parcorrf)>0.0000001 & abs(old_parcorrf)~=1);
spareM=sparse(spareM);
genenames=strrep(cellstr( deblank(casenames) ),'"','');
genenames(2)={'ste3'};
genenames(3)={'ste4'};
genenames(6)={'bni1'};
genenames(10)={'fus3'};
genenames(12)={'rst1'};
genenames(13)={'rst2'};
[x,y,h]=draw_graph(spareM,lower(genenames));
%test HOG path way %qulte a lot of missing data
clear all
close all
filename='hog_pathway.dat';
[tdata, varnames, casenames] = tblread(filename,'\t') ;
%t_data=tdata([1:7,9:11],:);
t_data=tdata;
t_data=som_normalize(t_data','var');
N=46;
q=13;
cutoff=0.05
[old_parcorrf,delet_edge,record_edge]=MGraph_GaussModel(cov(t_data),N,q,cutoff)
spareM=(abs(old_parcorrf)>0.0000001 & abs(old_parcorrf)~=1);
spareM=sparse(spareM);
genenames=strrep(cellstr( deblank(casenames) ),'"','');
genenames(1)={'sln1'};
genenames(2)={'sho1'};
genenames(9)={'pbs2'};
genenames(10)={'hog1'};
genenames(13)={'mcm1'};
[x,y,h]=draw_graph(spareM,lower(genenames));
%test Filamentous path way %qulte a lot of missing data
%should include all the genes, there are some missing values
%something wrong with the draw graph, if there is a line cross the middle of the cirle
%that means this nodes had correlation with the nearest nodes with no line cross the circle
clear all
close all
filename='filamentous_pathway.dat';
[tdata, varnames, casenames] = tblread(filename,'\t') ;
%t_data=tdata([1:7,9,11],:);
t_data=tdata;
t_data=som_normalize(t_data','var');
N=46;
q=11;
cutoff=0.05
[old_parcorrf,delet_edge,record_edge]=MGraph_GaussModel(cov(t_data),N,q,cutoff)
spareM=(abs(old_parcorrf)>0.0000001 & abs(old_parcorrf)~=1);
spareM=sparse(spareM);
genenames=strrep(cellstr( deblank(casenames) ),'"','');
genenames(1)={'Sho1'};
genenames(2)={'Ras2'};
genenames(8)={'Rst1'};
genenames(9)={'Rst2'};
genenames(11)={'Tec1'};
[x,y,h]=draw_graph(spareM,lower(genenames));
