24 private final head_t[] head;
27 Cache(
int l_,
long size_)
32 for(
int i=0;i<l;i++) head[i] =
new head_t();
35 size = Math.max(size, 2* (
long) l);
37 lru_head.next = lru_head.prev = lru_head;
40 private void lru_delete(head_t h)
47 private void lru_insert(head_t h)
51 h.prev = lru_head.prev;
60 int get_data(
int index,
float[][] data,
int len)
62 head_t h = head[
index];
63 if(h.len > 0) lru_delete(h);
64 int more = len - h.len;
71 head_t old = lru_head.next;
79 float[] new_data =
new float[len];
80 if(h.data !=
null) System.arraycopy(h.data,0,new_data,0,h.len);
83 do {
int _=h.len; h.len=len; len=_;}
while(
false);
91 void swap_index(
int i,
int j)
95 if(head[i].len > 0) lru_delete(head[i]);
96 if(head[j].len > 0) lru_delete(head[j]);
97 do {
float[] _=head[i].data; head[i].data=head[j].data; head[j].data=_;}
while(
false);
98 do {
int _=head[i].len; head[i].len=head[j].len; head[j].len=_;}
while(
false);
99 if(head[i].len > 0) lru_insert(head[i]);
100 if(head[j].len > 0) lru_insert(head[j]);
102 if(i>j)
do {
int _=i; i=j; j=_;}
while(
false);
103 for(head_t h = lru_head.next; h!=lru_head; h=h.next)
108 do {
float _=h.data[i]; h.data[i]=h.data[j]; h.data[j]=_;}
while(
false);
130 abstract float[] get_Q(
int column,
int len);
131 abstract double[] get_QD();
132 abstract void swap_index(
int i,
int j);
137 private final double[] x_square;
140 private final int kernel_type;
141 private final int degree;
142 private final double gamma;
143 private final double coef0;
145 abstract float[] get_Q(
int column,
int len);
146 abstract double[] get_QD();
148 void swap_index(
int i,
int j)
150 do {
svm_node[] _=
x[i];
x[i]=
x[j];
x[j]=_;}
while(
false);
151 if(x_square !=
null)
do {
double _=x_square[i]; x_square[i]=x_square[j]; x_square[j]=_;}
while(
false);
154 private static double powi(
double base,
int times)
156 double tmp = base, ret = 1.0;
158 for(
int t=times; t>0; t/=2)
166 double kernel_function(
int i,
int j)
171 return dot(
x[i],
x[j]);
173 return powi(gamma*
dot(
x[i],
x[j])+coef0,degree);
175 return Math.exp(-gamma*(x_square[i]+x_square[j]-2*
dot(
x[i],
x[j])));
177 return Math.tanh(gamma*
dot(
x[i],
x[j])+coef0);
179 return x[i][(int)(
x[j][0].value)].
value;
196 x_square =
new double[l];
198 x_square[i] =
dot(
x[i],
x[i]);
200 else x_square =
null;
210 while(i < xlen && j < ylen)
241 while(i < xlen && j < ylen)
277 return x[(int)(y[0].value)].
value;
306 static final byte LOWER_BOUND = 0;
307 static final byte UPPER_BOUND = 1;
308 static final byte FREE = 2;
321 static final double INF = java.lang.Double.POSITIVE_INFINITY;
325 return (y[i] > 0)? Cp : Cn;
327 void update_alpha_status(
int i)
329 if(alpha[i] >= get_C(i))
330 alpha_status[i] = UPPER_BOUND;
331 else if(alpha[i] <= 0)
332 alpha_status[i] = LOWER_BOUND;
333 else alpha_status[i] = FREE;
335 boolean is_upper_bound(
int i) {
return alpha_status[i] == UPPER_BOUND; }
336 boolean is_lower_bound(
int i) {
return alpha_status[i] == LOWER_BOUND; }
337 boolean is_free(
int i) {
return alpha_status[i] == FREE; }
341 static class SolutionInfo {
344 double upper_bound_p;
345 double upper_bound_n;
349 void swap_index(
int i,
int j)
352 do {
byte _=y[i]; y[i]=y[j]; y[j]=_;}
while(
false);
353 do {
double _=G[i]; G[i]=G[j]; G[j]=_;}
while(
false);
354 do {
byte _=alpha_status[i]; alpha_status[i]=alpha_status[j]; alpha_status[j]=_;}
while(
false);
355 do {
double _=alpha[i]; alpha[i]=alpha[j]; alpha[j]=_;}
while(
false);
356 do {
double _=p[i]; p[i]=p[j]; p[j]=_;}
while(
false);
357 do {
int _=active_set[i]; active_set[i]=active_set[j]; active_set[j]=_;}
while(
false);
358 do {
double _=G_bar[i]; G_bar[i]=G_bar[j]; G_bar[j]=_;}
while(
false);
361 void reconstruct_gradient()
365 if(active_size == l)
return;
370 for(j=active_size;j<l;j++)
371 G[j] = G_bar[j] + p[j];
373 for(j=0;j<active_size;j++)
377 if(2*nr_free < active_size)
378 svm.info(
"\nWARNING: using -h 0 may be faster\n");
380 if (nr_free*l > 2*active_size*(l-active_size))
382 for(i=active_size;i<l;i++)
384 float[] Q_i = Q.
get_Q(i,active_size);
385 for(j=0;j<active_size;j++)
387 G[i] += alpha[j] * Q_i[j];
392 for(i=0;i<active_size;i++)
395 float[] Q_i = Q.
get_Q(i,l);
396 double alpha_i = alpha[i];
397 for(j=active_size;j<l;j++)
398 G[j] += alpha_i * Q_i[j];
403 void Solve(
int l,
QMatrix Q,
double[] p_,
byte[] y_,
404 double[] alpha_,
double Cp,
double Cn,
double eps, SolutionInfo si,
int shrinking)
409 p = (
double[])p_.clone();
410 y = (
byte[])y_.clone();
411 alpha = (
double[])alpha_.clone();
415 this.unshrink =
false;
419 alpha_status =
new byte[l];
421 update_alpha_status(i);
426 active_set =
new int[l];
435 G_bar =
new double[l];
443 if(!is_lower_bound(i))
445 float[] Q_i = Q.
get_Q(i,l);
446 double alpha_i = alpha[i];
449 G[j] += alpha_i*Q_i[j];
450 if(is_upper_bound(i))
452 G_bar[j] += get_C(i) * Q_i[j];
459 int max_iter = Math.max(10000000, l>Integer.MAX_VALUE/100 ? Integer.MAX_VALUE : 100*l);
460 int counter = Math.min(l,1000)+1;
461 int[] working_set =
new int[2];
463 while(iter < max_iter)
469 counter = Math.min(l,1000);
470 if(shrinking!=0) do_shrinking();
474 if(select_working_set(working_set)!=0)
477 reconstruct_gradient();
481 if(select_working_set(working_set)!=0)
487 int i = working_set[0];
488 int j = working_set[1];
494 float[] Q_i = Q.
get_Q(i,active_size);
495 float[] Q_j = Q.
get_Q(j,active_size);
497 double C_i = get_C(i);
498 double C_j = get_C(j);
500 double old_alpha_i = alpha[i];
501 double old_alpha_j = alpha[j];
505 double quad_coef = QD[i]+QD[j]+2*Q_i[j];
508 double delta = (-G[i]-G[j])/quad_coef;
509 double diff = alpha[i] - alpha[j];
534 alpha[j] = C_i - diff;
542 alpha[i] = C_j + diff;
548 double quad_coef = QD[i]+QD[j]-2*Q_i[j];
551 double delta = (G[i]-G[j])/quad_coef;
552 double sum = alpha[i] + alpha[j];
561 alpha[j] = sum - C_i;
577 alpha[i] = sum - C_j;
592 double delta_alpha_i = alpha[i] - old_alpha_i;
593 double delta_alpha_j = alpha[j] - old_alpha_j;
595 for(
int k=0;k<active_size;k++)
597 G[k] += Q_i[k]*delta_alpha_i + Q_j[k]*delta_alpha_j;
603 boolean ui = is_upper_bound(i);
604 boolean uj = is_upper_bound(j);
605 update_alpha_status(i);
606 update_alpha_status(j);
608 if(ui != is_upper_bound(i))
613 G_bar[k] -= C_i * Q_i[k];
616 G_bar[k] += C_i * Q_i[k];
619 if(uj != is_upper_bound(j))
624 G_bar[k] -= C_j * Q_j[k];
627 G_bar[k] += C_j * Q_j[k];
638 reconstruct_gradient();
642 svm.info(
"\nWARNING: reaching max number of iterations");
647 si.rho = calculate_rho();
654 v += alpha[i] * (G[i] + p[i]);
662 alpha_[active_set[i]] = alpha[i];
665 si.upper_bound_p = Cp;
666 si.upper_bound_n = Cn;
668 svm.info(
"\noptimization finished, #iter = "+iter+
"\n");
672 int select_working_set(
int[] working_set)
684 double obj_diff_min =
INF;
686 for(
int t=0;t<active_size;t++)
689 if(!is_upper_bound(t))
698 if(!is_lower_bound(t))
709 Q_i = Q.
get_Q(i,active_size);
711 for(
int j=0;j<active_size;j++)
715 if (!is_lower_bound(j))
717 double grad_diff=Gmax+G[j];
723 double quad_coef = QD[i]+QD[j]-2.0*y[i]*Q_i[j];
725 obj_diff = -(grad_diff*grad_diff)/quad_coef;
727 obj_diff = -(grad_diff*grad_diff)/1e-12;
729 if (obj_diff <= obj_diff_min)
732 obj_diff_min = obj_diff;
739 if (!is_upper_bound(j))
741 double grad_diff= Gmax-G[j];
747 double quad_coef = QD[i]+QD[j]+2.0*y[i]*Q_i[j];
749 obj_diff = -(grad_diff*grad_diff)/quad_coef;
751 obj_diff = -(grad_diff*grad_diff)/1e-12;
753 if (obj_diff <= obj_diff_min)
756 obj_diff_min = obj_diff;
766 working_set[0] = Gmax_idx;
767 working_set[1] = Gmin_idx;
771 private boolean be_shrunk(
int i,
double Gmax1,
double Gmax2)
773 if(is_upper_bound(i))
776 return(-G[i] > Gmax1);
778 return(-G[i] > Gmax2);
780 else if(is_lower_bound(i))
783 return(G[i] > Gmax2);
785 return(G[i] > Gmax1);
798 for(i=0;i<active_size;i++)
802 if(!is_upper_bound(i))
807 if(!is_lower_bound(i))
815 if(!is_upper_bound(i))
820 if(!is_lower_bound(i))
828 if(unshrink ==
false && Gmax1 + Gmax2 <= eps*10)
831 reconstruct_gradient();
835 for(i=0;i<active_size;i++)
836 if (be_shrunk(i, Gmax1, Gmax2))
839 while (active_size > i)
841 if (!be_shrunk(active_size, Gmax1, Gmax2))
843 swap_index(i,active_size);
851 double calculate_rho()
855 double ub =
INF, lb = -
INF, sum_free = 0;
856 for(
int i=0;i<active_size;i++)
858 double yG = y[i]*G[i];
860 if(is_lower_bound(i))
863 ub = Math.min(ub,yG);
865 lb = Math.max(lb,yG);
867 else if(is_upper_bound(i))
870 ub = Math.min(ub,yG);
872 lb = Math.max(lb,yG);
882 r = sum_free/nr_free;
898 private SolutionInfo si;
900 void Solve(
int l,
QMatrix Q,
double[] p,
byte[] y,
901 double[] alpha,
double Cp,
double Cn,
double eps,
902 SolutionInfo si,
int shrinking)
905 super.Solve(l,Q,p,y,alpha,Cp,Cn,eps,si,shrinking);
909 int select_working_set(
int[] working_set)
918 double Gmaxp2 = -
INF;
922 double Gmaxn2 = -
INF;
926 double obj_diff_min =
INF;
928 for(
int t=0;t<active_size;t++)
931 if(!is_upper_bound(t))
940 if(!is_lower_bound(t))
953 Q_ip = Q.
get_Q(ip,active_size);
955 Q_in = Q.
get_Q(in,active_size);
957 for(
int j=0;j<active_size;j++)
961 if (!is_lower_bound(j))
963 double grad_diff=Gmaxp+G[j];
969 double quad_coef = QD[ip]+QD[j]-2*Q_ip[j];
971 obj_diff = -(grad_diff*grad_diff)/quad_coef;
973 obj_diff = -(grad_diff*grad_diff)/1e-12;
975 if (obj_diff <= obj_diff_min)
978 obj_diff_min = obj_diff;
985 if (!is_upper_bound(j))
987 double grad_diff=Gmaxn-G[j];
993 double quad_coef = QD[in]+QD[j]-2*Q_in[j];
995 obj_diff = -(grad_diff*grad_diff)/quad_coef;
997 obj_diff = -(grad_diff*grad_diff)/1e-12;
999 if (obj_diff <= obj_diff_min)
1002 obj_diff_min = obj_diff;
1009 if(Math.max(Gmaxp+Gmaxp2,Gmaxn+Gmaxn2) < eps)
1012 if(y[Gmin_idx] == +1)
1013 working_set[0] = Gmaxp_idx;
1015 working_set[0] = Gmaxn_idx;
1016 working_set[1] = Gmin_idx;
1021 private boolean be_shrunk(
int i,
double Gmax1,
double Gmax2,
double Gmax3,
double Gmax4)
1023 if(is_upper_bound(i))
1026 return(-G[i] > Gmax1);
1028 return(-G[i] > Gmax4);
1030 else if(is_lower_bound(i))
1033 return(G[i] > Gmax2);
1035 return(G[i] > Gmax3);
1043 double Gmax1 = -
INF;
1044 double Gmax2 = -
INF;
1045 double Gmax3 = -
INF;
1046 double Gmax4 = -
INF;
1050 for(i=0;i<active_size;i++)
1052 if(!is_upper_bound(i))
1056 if(-G[i] > Gmax1) Gmax1 = -G[i];
1058 else if(-G[i] > Gmax4) Gmax4 = -G[i];
1060 if(!is_lower_bound(i))
1064 if(G[i] > Gmax2) Gmax2 = G[i];
1066 else if(G[i] > Gmax3) Gmax3 = G[i];
1070 if(unshrink ==
false && Math.max(Gmax1+Gmax2,Gmax3+Gmax4) <= eps*10)
1073 reconstruct_gradient();
1077 for(i=0;i<active_size;i++)
1078 if (be_shrunk(i, Gmax1, Gmax2, Gmax3, Gmax4))
1081 while (active_size > i)
1083 if (!be_shrunk(active_size, Gmax1, Gmax2, Gmax3, Gmax4))
1085 swap_index(i,active_size);
1093 double calculate_rho()
1095 int nr_free1 = 0,nr_free2 = 0;
1096 double ub1 =
INF, ub2 =
INF;
1097 double lb1 = -
INF, lb2 = -
INF;
1098 double sum_free1 = 0, sum_free2 = 0;
1100 for(
int i=0;i<active_size;i++)
1104 if(is_lower_bound(i))
1105 ub1 = Math.min(ub1,G[i]);
1106 else if(is_upper_bound(i))
1107 lb1 = Math.max(lb1,G[i]);
1116 if(is_lower_bound(i))
1117 ub2 = Math.min(ub2,G[i]);
1118 else if(is_upper_bound(i))
1119 lb2 = Math.max(lb2,G[i]);
1130 r1 = sum_free1/nr_free1;
1135 r2 = sum_free2/nr_free2;
1149 private final byte[] y;
1150 private final Cache cache;
1151 private final double[] QD;
1156 y = (
byte[])y_.clone();
1158 QD =
new double[
prob.
l];
1159 for(
int i=0;i<
prob.
l;i++)
1160 QD[i] = kernel_function(i,i);
1163 float[] get_Q(
int i,
int len)
1165 float[][] data =
new float[1][];
1167 if((start = cache.get_data(i,data,len)) < len)
1169 for(j=start;j<len;j++)
1170 data[0][j] = (
float)(y[i]*y[j]*kernel_function(i,j));
1180 void swap_index(
int i,
int j)
1182 cache.swap_index(i,j);
1183 super.swap_index(i,j);
1184 do {
byte _=y[i]; y[i]=y[j]; y[j]=_;}
while(
false);
1185 do {
double _=QD[i]; QD[i]=QD[j]; QD[j]=_;}
while(
false);
1191 private final Cache cache;
1192 private final double[] QD;
1198 QD =
new double[
prob.
l];
1199 for(
int i=0;i<
prob.
l;i++)
1200 QD[i] = kernel_function(i,i);
1203 float[] get_Q(
int i,
int len)
1205 float[][] data =
new float[1][];
1207 if((start = cache.get_data(i,data,len)) < len)
1209 for(j=start;j<len;j++)
1210 data[0][j] = (
float)kernel_function(i,j);
1220 void swap_index(
int i,
int j)
1222 cache.swap_index(i,j);
1223 super.swap_index(i,j);
1224 do {
double _=QD[i]; QD[i]=QD[j]; QD[j]=_;}
while(
false);
1230 private final int l;
1231 private final Cache cache;
1232 private final byte[] sign;
1233 private final int[]
index;
1234 private int next_buffer;
1235 private float[][]
buffer;
1236 private final double[] QD;
1243 QD =
new double[2*l];
1244 sign =
new byte[2*l];
1245 index =
new int[2*l];
1246 for(
int k=0;k<l;k++)
1252 QD[k] = kernel_function(k,k);
1255 buffer =
new float[2][2*l];
1259 void swap_index(
int i,
int j)
1261 do {
byte _=sign[i]; sign[i]=sign[j]; sign[j]=_;}
while(
false);
1263 do {
double _=QD[i]; QD[i]=QD[j]; QD[j]=_;}
while(
false);
1266 float[] get_Q(
int i,
int len)
1268 float[][] data =
new float[1][];
1269 int j, real_i =
index[i];
1270 if(cache.get_data(real_i,data,l) < l)
1273 data[0][j] = (
float)kernel_function(real_i,j);
1277 float buf[] =
buffer[next_buffer];
1278 next_buffer = 1 - next_buffer;
1281 buf[j] = (
float) si * sign[j] * data[0][
index[j]];
1296 public static final Random
rand =
new Random();
1300 public void print(String
s)
1302 System.out.print(
s);
1309 static void info(String s)
1315 double[] alpha, Solver.SolutionInfo si,
1316 double Cp,
double Cn)
1319 double[] minus_ones =
new double[l];
1320 byte[] y =
new byte[l];
1328 if(
prob.
y[i] > 0) y[i] = +1;
else y[i] = -1;
1331 Solver
s =
new Solver();
1332 s.Solve(l,
new SVC_Q(
prob,
param,y), minus_ones, y,
1337 sum_alpha += alpha[i];
1340 svm.info(
"nu = "+sum_alpha/(Cp*
prob.
l)+
"\n");
1347 double[] alpha, Solver.SolutionInfo si)
1353 byte[] y =
new byte[l];
1361 double sum_pos = nu*l/2;
1362 double sum_neg = nu*l/2;
1367 alpha[i] = Math.min(1.0,sum_pos);
1368 sum_pos -= alpha[i];
1372 alpha[i] = Math.min(1.0,sum_neg);
1373 sum_neg -= alpha[i];
1376 double[] zeros =
new double[l];
1381 Solver_NU
s =
new Solver_NU();
1382 s.Solve(l,
new SVC_Q(
prob,
param,y), zeros, y,
1386 svm.info(
"C = "+1/r+
"\n");
1393 si.upper_bound_p = 1/r;
1394 si.upper_bound_n = 1/r;
1398 double[] alpha, Solver.SolutionInfo si)
1401 double[] zeros =
new double[l];
1402 byte[] ones =
new byte[l];
1420 Solver
s =
new Solver();
1421 s.Solve(l,
new ONE_CLASS_Q(
prob,
param), zeros, ones,
1426 double[] alpha, Solver.SolutionInfo si)
1429 double[] alpha2 =
new double[2*l];
1430 double[] linear_term =
new double[2*l];
1431 byte[] y =
new byte[2*l];
1445 Solver
s =
new Solver();
1446 s.Solve(2*l,
new SVR_Q(
prob,
param), linear_term, y,
1449 double sum_alpha = 0;
1452 alpha[i] = alpha2[i] - alpha2[i+l];
1453 sum_alpha += Math.abs(alpha[i]);
1455 svm.info(
"nu = "+sum_alpha/(
param.
C*l)+
"\n");
1459 double[] alpha, Solver.SolutionInfo si)
1463 double[] alpha2 =
new double[2*l];
1464 double[] linear_term =
new double[2*l];
1465 byte[] y =
new byte[2*l];
1468 double sum = C *
param.
nu * l / 2;
1471 alpha2[i] = alpha2[i+l] = Math.min(sum,C);
1474 linear_term[i] = -
prob.
y[i];
1477 linear_term[i+l] =
prob.
y[i];
1481 Solver_NU
s =
new Solver_NU();
1482 s.Solve(2*l,
new SVR_Q(
prob,
param), linear_term, y,
1485 svm.info(
"epsilon = "+(-si.r)+
"\n");
1488 alpha[i] = alpha2[i] - alpha2[i+l];
1502 double Cp,
double Cn)
1504 double[] alpha =
new double[
prob.
l];
1525 svm.info(
"obj = "+si.obj+
", rho = "+si.rho+
"\n");
1531 for(
int i=0;i<
prob.
l;i++)
1533 if(Math.abs(alpha[i]) > 0)
1538 if(Math.abs(alpha[i]) >= si.upper_bound_p)
1543 if(Math.abs(alpha[i]) >= si.upper_bound_n)
1549 svm.info(
"nSV = "+nSV+
", nBSV = "+nBSV+
"\n");
1562 double prior1=0, prior0 = 0;
1566 if (labels[i] > 0) prior1+=1;
1570 double min_step=1e-10;
1573 double hiTarget=(prior1+1.0)/(prior1+2.0);
1574 double loTarget=1/(prior0+2.0);
1575 double[] t=
new double[l];
1576 double fApB,p,q,h11,h22,h21,g1,g2,det,dA,dB,gd,stepsize;
1577 double newA,newB,newf,d1,d2;
1581 A=0.0; B=Math.log((prior0+1.0)/(prior1+1.0));
1586 if (labels[i]>0) t[i]=hiTarget;
1588 fApB = dec_values[i]*A+B;
1590 fval += t[i]*fApB + Math.log(1+Math.exp(-fApB));
1592 fval += (t[i] - 1)*fApB +Math.log(1+Math.exp(fApB));
1594 for (iter=0;iter<max_iter;iter++)
1599 h21=0.0;g1=0.0;g2=0.0;
1602 fApB = dec_values[i]*A+B;
1605 p=Math.exp(-fApB)/(1.0+Math.exp(-fApB));
1606 q=1.0/(1.0+Math.exp(-fApB));
1610 p=1.0/(1.0+Math.exp(fApB));
1611 q=Math.exp(fApB)/(1.0+Math.exp(fApB));
1614 h11+=dec_values[i]*dec_values[i]*d2;
1616 h21+=dec_values[i]*d2;
1618 g1+=dec_values[i]*d1;
1623 if (Math.abs(g1)<eps && Math.abs(g2)<eps)
1627 det=h11*h22-h21*h21;
1628 dA=-(h22*g1 - h21 * g2) / det;
1629 dB=-(-h21*g1+ h11 * g2) / det;
1634 while (stepsize >= min_step)
1636 newA = A + stepsize * dA;
1637 newB = B + stepsize * dB;
1643 fApB = dec_values[i]*newA+newB;
1645 newf += t[i]*fApB + Math.log(1+Math.exp(-fApB));
1647 newf += (t[i] - 1)*fApB +Math.log(1+Math.exp(fApB));
1650 if (newf<fval+0.0001*stepsize*gd)
1652 A=newA;B=newB;fval=newf;
1656 stepsize = stepsize / 2.0;
1659 if (stepsize < min_step)
1661 svm.info(
"Line search fails in two-class probability estimates\n");
1667 svm.info(
"Reaching maximal iterations in two-class probability estimates\n");
1668 probAB[0]=A;probAB[1]=B;
1673 double fApB = decision_value*A+B;
1675 return Math.exp(-fApB)/(1.0+Math.exp(-fApB));
1677 return 1.0/(1+Math.exp(fApB)) ;
1684 int iter = 0, max_iter=Math.max(100,k);
1685 double[][] Q=
new double[k][k];
1686 double[] Qp=
new double[k];
1687 double pQp, eps=0.005/k;
1695 Q[t][t]+=r[j][t]*r[j][t];
1700 Q[t][t]+=r[j][t]*r[j][t];
1701 Q[t][j]=-r[j][t]*r[t][j];
1704 for (iter=0;iter<max_iter;iter++)
1712 Qp[t]+=Q[t][j]*p[j];
1718 double error=Math.abs(Qp[t]-pQp);
1719 if (error>max_error)
1722 if (max_error<eps)
break;
1726 double diff=(-Qp[t]+pQp)/Q[t][t];
1728 pQp=(pQp+diff*(diff*Q[t][t]+2*Qp[t]))/(1+diff)/(1+diff);
1731 Qp[j]=(Qp[j]+diff*Q[t][j])/(1+diff);
1737 svm.info(
"Exceeds max_iter in multiclass_prob\n");
1745 int[] perm =
new int[
prob.
l];
1746 double[] dec_values =
new double[
prob.
l];
1749 for(i=0;i<
prob.
l;i++) perm[i]=i;
1750 for(i=0;i<
prob.
l;i++)
1753 do {
int _=perm[i]; perm[i]=perm[j]; perm[j]=_;}
while(
false);
1762 subprob.
l =
prob.
l-(end-begin);
1763 subprob.x =
new svm_node[subprob.l][];
1764 subprob.y =
new double[subprob.l];
1767 for(j=0;j<begin;j++)
1769 subprob.x[k] =
prob.
x[perm[j]];
1770 subprob.y[k] =
prob.
y[perm[j]];
1773 for(j=end;j<
prob.
l;j++)
1775 subprob.x[k] =
prob.
x[perm[j]];
1776 subprob.y[k] =
prob.
y[perm[j]];
1779 int p_count=0,n_count=0;
1786 if(p_count==0 && n_count==0)
1787 for(j=begin;j<end;j++)
1788 dec_values[perm[j]] = 0;
1789 else if(p_count > 0 && n_count == 0)
1790 for(j=begin;j<end;j++)
1791 dec_values[perm[j]] = 1;
1792 else if(p_count == 0 && n_count > 0)
1793 for(j=begin;j<end;j++)
1794 dec_values[perm[j]] = -1;
1802 subparam.
weight =
new double[2];
1808 for(j=begin;j<end;j++)
1810 double[] dec_value=
new double[1];
1812 dec_values[perm[j]]=dec_value[0];
1814 dec_values[perm[j]] *= submodel.
label[0];
1826 double[] ymv =
new double[
prob.
l];
1832 for(i=0;i<
prob.
l;i++)
1834 ymv[i]=
prob.
y[i]-ymv[i];
1835 mae += Math.abs(ymv[i]);
1838 double std=Math.sqrt(2*mae*mae);
1841 for(i=0;i<
prob.
l;i++)
1842 if (Math.abs(ymv[i]) > 5*
std)
1845 mae+=Math.abs(ymv[i]);
1846 mae /= (
prob.
l-count);
1847 svm.info(
"Prob. model for test data: target value = predicted value + z,\nz: Laplace distribution e^(-|z|/sigma)/(2sigma),sigma="+mae+
"\n");
1856 int max_nr_class = 16;
1858 int[]
label =
new int[max_nr_class];
1859 int[] count =
new int[max_nr_class];
1860 int[] data_label =
new int[l];
1865 int this_label = (int)(
prob.
y[i]);
1867 for(j=0;j<nr_class;j++)
1869 if(this_label ==
label[j])
1878 if(nr_class == max_nr_class)
1881 int[] new_data =
new int[max_nr_class];
1882 System.arraycopy(
label,0,new_data,0,
label.length);
1884 new_data =
new int[max_nr_class];
1885 System.arraycopy(count,0,new_data,0,count.length);
1888 label[nr_class] = this_label;
1889 count[nr_class] = 1;
1894 int[]
start =
new int[nr_class];
1896 for(i=1;i<nr_class;i++)
1900 perm[
start[data_label[i]]] = i;
1901 ++
start[data_label[i]];
1904 for(i=1;i<nr_class;i++)
1907 nr_class_ret[0] = nr_class;
1908 label_ret[0] =
label;
1909 start_ret[0] =
start;
1910 count_ret[0] = count;
1940 decision_function
f = svm_train_one(
prob,
param,0,0);
1946 for(i=0;i<
prob.
l;i++)
1947 if(Math.abs(
f.alpha[i]) > 0) ++nSV;
1952 for(i=0;i<
prob.
l;i++)
1953 if(Math.abs(
f.alpha[i]) > 0)
1964 int[] tmp_nr_class =
new int[1];
1965 int[][] tmp_label =
new int[1][];
1966 int[][] tmp_start =
new int[1][];
1967 int[][] tmp_count =
new int[1][];
1968 int[] perm =
new int[l];
1972 int nr_class = tmp_nr_class[0];
1973 int[]
label = tmp_label[0];
1974 int[]
start = tmp_start[0];
1975 int[] count = tmp_count[0];
1978 svm.info(
"WARNING: training data in only one class. See README for details.\n");
1987 double[] weighted_C =
new double[nr_class];
1988 for(i=0;i<nr_class;i++)
1993 for(j=0;j<nr_class;j++)
1997 System.err.print(
"WARNING: class label "+
param.
weight_label[i]+
" specified in weight is not found\n");
2004 boolean[] nonzero =
new boolean[l];
2007 decision_function[]
f =
new decision_function[nr_class*(nr_class-1)/2];
2009 double[] probA=
null,probB=
null;
2012 probA=
new double[nr_class*(nr_class-1)/2];
2013 probB=
new double[nr_class*(nr_class-1)/2];
2017 for(i=0;i<nr_class;i++)
2018 for(
int j=i+1;j<nr_class;j++)
2022 int ci = count[i], cj = count[j];
2025 sub_prob.
y =
new double[sub_prob.
l];
2029 sub_prob.
x[k] =
x[si+k];
2034 sub_prob.
x[ci+k] =
x[sj+k];
2035 sub_prob.
y[ci+k] = -1;
2040 double[] probAB=
new double[2];
2046 f[p] = svm_train_one(sub_prob,
param,weighted_C[i],weighted_C[j]);
2048 if(!nonzero[si+k] && Math.abs(
f[p].alpha[k]) > 0)
2049 nonzero[si+k] =
true;
2051 if(!nonzero[sj+k] && Math.abs(
f[p].alpha[ci+k]) > 0)
2052 nonzero[sj+k] =
true;
2061 for(i=0;i<nr_class;i++)
2064 model.
rho =
new double[nr_class*(nr_class-1)/2];
2065 for(i=0;i<nr_class*(nr_class-1)/2;i++)
2070 model.
probA =
new double[nr_class*(nr_class-1)/2];
2071 model.
probB =
new double[nr_class*(nr_class-1)/2];
2072 for(i=0;i<nr_class*(nr_class-1)/2;i++)
2085 int[] nz_count =
new int[nr_class];
2087 for(i=0;i<nr_class;i++)
2090 for(
int j=0;j<count[i];j++)
2091 if(nonzero[
start[i]+j])
2100 svm.info(
"Total nSV = "+nnz+
"\n");
2106 if(nonzero[i])
model.
SV[p++] =
x[i];
2108 int[] nz_start =
new int[nr_class];
2110 for(i=1;i<nr_class;i++)
2111 nz_start[i] = nz_start[i-1]+nz_count[i-1];
2114 for(i=0;i<nr_class-1;i++)
2118 for(i=0;i<nr_class;i++)
2119 for(
int j=i+1;j<nr_class;j++)
2130 int q = nz_start[i];
2149 int[] fold_start =
new int[
nr_fold+1];
2151 int[] perm =
new int[l];
2158 int[] tmp_nr_class =
new int[1];
2159 int[][] tmp_label =
new int[1][];
2160 int[][] tmp_start =
new int[1][];
2161 int[][] tmp_count =
new int[1][];
2165 int nr_class = tmp_nr_class[0];
2166 int[]
start = tmp_start[0];
2167 int[] count = tmp_count[0];
2170 int[] fold_count =
new int[
nr_fold];
2172 int[]
index =
new int[l];
2175 for (
c=0;
c<nr_class;
c++)
2176 for(i=0;i<count[
c];i++)
2178 int j = i+
rand.nextInt(count[
c]-i);
2184 for (
c=0;
c<nr_class;
c++)
2189 fold_start[i] = fold_start[i-1]+fold_count[i-1];
2190 for (
c=0;
c<nr_class;
c++)
2195 for(
int j=begin;j<end;j++)
2197 perm[fold_start[i]] =
index[j];
2203 fold_start[i] = fold_start[i-1]+fold_count[i-1];
2207 for(i=0;i<l;i++) perm[i]=i;
2210 int j = i+
rand.nextInt(l-i);
2211 do {
int _=perm[i]; perm[i]=perm[j]; perm[j]=_;}
while(
false);
2219 int begin = fold_start[i];
2220 int end = fold_start[i+1];
2224 subprob.
l = l-(end-begin);
2226 subprob.
y =
new double[subprob.
l];
2229 for(j=0;j<begin;j++)
2231 subprob.
x[k] =
prob.
x[perm[j]];
2232 subprob.
y[k] =
prob.
y[perm[j]];
2237 subprob.
x[k] =
prob.
x[perm[j]];
2238 subprob.
y[k] =
prob.
y[perm[j]];
2247 for(j=begin;j<end;j++)
2251 for(j=begin;j<end;j++)
2280 System.err.print(
"Model doesn't contain information for SVR probability inference\n");
2297 dec_values[0] = sum;
2300 return (sum>0)?1:-1;
2309 double[] kvalue =
new double[l];
2313 int[]
start =
new int[nr_class];
2315 for(i=1;i<nr_class;i++)
2318 int[] vote =
new int[nr_class];
2319 for(i=0;i<nr_class;i++)
2323 for(i=0;i<nr_class;i++)
2324 for(
int j=i+1;j<nr_class;j++)
2336 sum += coef1[si+k] * kvalue[si+k];
2338 sum += coef2[sj+k] * kvalue[sj+k];
2340 dec_values[p] = sum;
2342 if(dec_values[p] > 0)
2349 int vote_max_idx = 0;
2350 for(i=1;i<nr_class;i++)
2351 if(vote[i] > vote[vote_max_idx])
2361 double[] dec_values;
2365 dec_values =
new double[1];
2367 dec_values =
new double[nr_class*(nr_class-1)/2];
2379 double[] dec_values =
new double[nr_class*(nr_class-1)/2];
2382 double min_prob=1e-7;
2383 double[][] pairwise_prob=
new double[nr_class][nr_class];
2386 for(i=0;i<nr_class;i++)
2387 for(
int j=i+1;j<nr_class;j++)
2390 pairwise_prob[j][i]=1-pairwise_prob[i][j];
2395 int prob_max_idx = 0;
2396 for(i=1;i<nr_class;i++)
2397 if(prob_estimates[i] > prob_estimates[prob_max_idx])
2405 static final String svm_type_table[] =
2407 "c_svc",
"nu_svc",
"one_class",
"epsilon_svr",
"nu_svr",
2410 static final String kernel_type_table[]=
2412 "linear",
"polynomial",
"rbf",
"sigmoid",
"precomputed"
2417 DataOutputStream fp =
new DataOutputStream(
new BufferedOutputStream(
new FileOutputStream(model_file_name)));
2421 fp.writeBytes(
"svm_type "+svm_type_table[
param.
svm_type]+
"\n");
2438 fp.writeBytes(
"nr_class "+nr_class+
"\n");
2439 fp.writeBytes(
"total_sv "+l+
"\n");
2442 fp.writeBytes(
"rho");
2443 for(
int i=0;i<nr_class*(nr_class-1)/2;i++)
2445 fp.writeBytes(
"\n");
2450 fp.writeBytes(
"label");
2451 for(
int i=0;i<nr_class;i++)
2453 fp.writeBytes(
"\n");
2458 fp.writeBytes(
"probA");
2459 for(
int i=0;i<nr_class*(nr_class-1)/2;i++)
2461 fp.writeBytes(
"\n");
2465 fp.writeBytes(
"probB");
2466 for(
int i=0;i<nr_class*(nr_class-1)/2;i++)
2468 fp.writeBytes(
"\n");
2473 fp.writeBytes(
"nr_sv");
2474 for(
int i=0;i<nr_class;i++)
2476 fp.writeBytes(
"\n");
2479 fp.writeBytes(
"SV\n");
2483 for(
int i=0;i<l;i++)
2485 for(
int j=0;j<nr_class-1;j++)
2486 fp.writeBytes(sv_coef[j][i]+
" ");
2490 fp.writeBytes(
"0:"+(
int)(p[0].
value));
2492 for(
int j=0;j<p.length;j++)
2493 fp.writeBytes(p[j].
index+
":"+p[j].
value+
" ");
2494 fp.writeBytes(
"\n");
2500 private static double atof(String s)
2502 return Double.valueOf(
s).doubleValue();
2507 return Integer.parseInt(
s);
2512 return svm_load_model(
new BufferedReader(
new FileReader(model_file_name)));
2530 String
cmd = fp.readLine();
2531 String arg =
cmd.substring(
cmd.indexOf(
' ')+1);
2533 if(
cmd.startsWith(
"svm_type"))
2536 for(i=0;i<svm_type_table.length;i++)
2538 if(arg.indexOf(svm_type_table[i])!=-1)
2544 if(i == svm_type_table.length)
2546 System.err.print(
"unknown svm type.\n");
2550 else if(
cmd.startsWith(
"kernel_type"))
2553 for(i=0;i<kernel_type_table.length;i++)
2555 if(arg.indexOf(kernel_type_table[i])!=-1)
2561 if(i == kernel_type_table.length)
2563 System.err.print(
"unknown kernel function.\n");
2567 else if(
cmd.startsWith(
"degree"))
2569 else if(
cmd.startsWith(
"gamma"))
2571 else if(
cmd.startsWith(
"coef0"))
2573 else if(
cmd.startsWith(
"nr_class"))
2575 else if(
cmd.startsWith(
"total_sv"))
2577 else if(
cmd.startsWith(
"rho"))
2581 StringTokenizer st =
new StringTokenizer(arg);
2582 for(
int i=0;i<n;i++)
2585 else if(
cmd.startsWith(
"label"))
2589 StringTokenizer st =
new StringTokenizer(arg);
2590 for(
int i=0;i<n;i++)
2593 else if(
cmd.startsWith(
"probA"))
2597 StringTokenizer st =
new StringTokenizer(arg);
2598 for(
int i=0;i<n;i++)
2601 else if(
cmd.startsWith(
"probB"))
2605 StringTokenizer st =
new StringTokenizer(arg);
2606 for(
int i=0;i<n;i++)
2609 else if(
cmd.startsWith(
"nr_sv"))
2613 StringTokenizer st =
new StringTokenizer(arg);
2614 for(
int i=0;i<n;i++)
2617 else if(
cmd.startsWith(
"SV"))
2623 System.err.print(
"unknown text in model file: ["+
cmd+
"]\n");
2635 for(
int i=0;i<l;i++)
2637 String
line = fp.readLine();
2638 StringTokenizer st =
new StringTokenizer(
line,
" \t\n\r\f:");
2640 for(
int k=0;k<m;k++)
2642 int n = st.countTokens()/2;
2644 for(
int j=0;j<n;j++)
2666 return "unknown svm type";
2676 return "unknown kernel type";
2682 return "degree of polynomial kernel < 0";
2687 return "cache_size <= 0";
2702 return "nu <= 0 or nu > 1";
2710 return "shrinking != 0 and shrinking != 1";
2714 return "probability != 0 and probability != 1";
2718 return "one-class SVM probability output not supported yet";
2725 int max_nr_class = 16;
2727 int[]
label =
new int[max_nr_class];
2728 int[] count =
new int[max_nr_class];
2733 int this_label = (int)
prob.
y[i];
2735 for(j=0;j<nr_class;j++)
2736 if(this_label ==
label[j])
2744 if(nr_class == max_nr_class)
2747 int[] new_data =
new int[max_nr_class];
2748 System.arraycopy(
label,0,new_data,0,
label.length);
2751 new_data =
new int[max_nr_class];
2752 System.arraycopy(count,0,new_data,0,count.length);
2755 label[nr_class] = this_label;
2756 count[nr_class] = 1;
2761 for(i=0;i<nr_class;i++)
2764 for(
int j=i+1;j<nr_class;j++)
2767 if(
param.
nu*(n1+n2)/2 > Math.min(n1,n2))
2768 return "specified nu is infeasible";
2789 if (print_func ==
null)