Training gets stuck on a specific dataset

[fumoboy007]

I am trying to train a C-SVC on a specific dataset. The training process gets stuck, never finishing.

Reproduction Steps

1. Unarchive libsvm-stuck-training-issue.zip, which includes reproduction code and data. (The data has already been scaled.)
2. Copy svm.h and svm.cpp from LIBSVM v335 (or any older version up to and including v300) to the same directory.
3. clang++ -std=c++11 -o reproduce_issue reproduce_issue.cpp svm.cpp
4. ./reproduce_issue

Preliminary Investigation

• The issue is not reproducible after reverting 1c80a42, though I’m not sure why. Similarly, changing typedef float Qfloat to typedef double Qfloat causes the issue to go away. This could be a red herring.
• After adding some logging to Solver::select_working_set, I can see that the training process eventually gets stuck on the working set {1525, 2023}, repeatedly selecting these two indices for i and j but with alternating order.

My hunch is that there is a bug in the working set selection or the stopping criteria. I will leave it to the experts to investigate further!

[cjlin1]

Your data is not a normal one as it has only 1 feature.
Using a high degree with a large regularization parameter easily
causes numerical difficulties.

I can reproduce your results:

$ make; ./svm-train -t 1 -e 0.01 -d 4 -c 100 -g 0.9597420397825849 -w0 1.04884106 -w1 0.95550528 ../test3/data1.txt
make: Nothing to be done for 'all'.
........................................................
WARNING: using -h 0 may be faster
*...............
WARNING: using -h 0 may be faster
*.................
WARNING: using -h 0 may be faster
*..........................................................................
WARNING: using -h 0 may be faster
*
WARNING: reaching max number of iterations

optimization finished, #iter = 10000000
obj = -244346.197253, rho = -1.203479
nSV = 2448, nBSV = 2440
Total nSV = 2448

However, if I change

typedef float Qfloat;

to

typedef double Qfloat;

then the problem disappears

$ make; ./svm-train -t 1 -e 0.01 -d 4 -c 100 -g 0.9597420397825849 -w0 1.04884106 -w1 0.95550528 ../test3/data1.txt
make: Nothing to be done for 'all'.
............................
WARNING: using -h 0 may be faster
*............................................
WARNING: using -h 0 may be faster
*.
WARNING: using -h 0 may be faster
*..................
WARNING: using -h 0 may be faster
*...........................
WARNING: using -h 0 may be faster
*.........
WARNING: using -h 0 may be faster
*............
WARNING: using -h 0 may be faster
*.............
WARNING: using -h 0 may be faster
*..............
WARNING: using -h 0 may be faster
*.....
WARNING: using -h 0 may be faster
*.................
WARNING: using -h 0 may be faster
*...........
WARNING: using -h 0 may be faster
*
optimization finished, #iter = 196403
obj = -245104.184532, rho = -1.003508
nSV = 2455, nBSV = 2448
Total nSV = 2455

If you diff 2.91 and 3.0, the main change is

$ diff ./svm.cpp ../libsvm-3.0/svm.cpp
195c195
< virtual Qfloat *get_QD() const = 0;

virtual double *get_QD() const = 0;
208c208
< virtual Qfloat *get_QD() const = 0;

---

virtual double *get_QD() const = 0;
415c415
< const Qfloat *QD;

---

const double QD;
601c601
< double quad_coef = Q_i[i]+Q_j[j]+2Q_i[j];

---

  	double quad_coef = QD[i]+QD[j]+2*Q_i[j];

644c644
< double quad_coef = Q_i[i]+Q_j[j]-2*Q_i[j];

  	double quad_coef = QD[i]+QD[j]-2*Q_i[j];

822c822
< double quad_coef=Q_i[i]+QD[j]-2.0*y[i]*Q_i[j];

  			double quad_coef = QD[i]+QD[j]-2.0*y[i]*Q_i[j];

846c846
< double quad_coef=Q_i[i]+QD[j]+2.0*y[i]*Q_i[j];

  			double quad_coef = QD[i]+QD[j]+2.0*y[i]*Q_i[j];

1074c1074
< double quad_coef = Q_ip[ip]+QD[j]-2*Q_ip[j];

  			double quad_coef = QD[ip]+QD[j]-2*Q_ip[j];

1098c1098
< double quad_coef = Q_in[in]+QD[j]-2*Q_in[j];

  			double quad_coef = QD[in]+QD[j]-2*Q_in[j];

1259c1259
< QD = new Qfloat[prob.l];

  QD = new double[prob.l];

1261c1261
< QD[i]= (Qfloat)(this->*kernel_function)(i,i);

  	QD[i] = (this->*kernel_function)(i,i);

1276c1276
< Qfloat *get_QD() const

double *get_QD() const
1298c1298
< Qfloat *QD;

---

double *QD;
1308c1308
< QD = new Qfloat[prob.l];

---

  QD = new double[prob.l];

1310c1310
< QD[i]= (Qfloat)(this->*kernel_function)(i,i);

  	QD[i] = (this->*kernel_function)(i,i);

1325c1325
< Qfloat *get_QD() const

double *get_QD() const
1344c1344
< Qfloat *QD;

---

double QD;
1355c1355
< QD = new Qfloat[2l];

---

  QD = new double[2*l];

1364,1365c1364,1365
< QD[k]= (Qfloat)(this->*kernel_function)(k,k);
< QD[k+l]=QD[k];

  	QD[k] = (this->*kernel_function)(k,k);
  	QD[k+l] = QD[k];

1398c1398
< Qfloat *get_QD() const

double *get_QD() const
1419c1419
< Qfloat *QD;

---

double *QD;
1439c1439

I believe for the case you tried, certain numerical issues occur.
According to our log, in 2.91 we did the change

QD is double now. All Q_i[i] are replaced by QD[i]
to avoid numerical problems.

It aims to solve some numerical problems.

In summary, numerical issues do occur in libsvm in different circumstances.
Yours happens to be one where the old setting works. But we did that
change because of some situations where the old one does not work.

[fumoboy007]

Hi @cjlin1, thanks for taking a look!

Is there a way that LIBSVM can detect numerical difficulties and respond accordingly? I mentioned in #225 (comment) that the working set selection keeps choosing the same two indices without making progress. This seems like a situation that can be detected and mitigated.

Your data is not a normal one as it has only 1 feature.

Using a high degree with a large regularization parameter easily causes numerical difficulties.

I can understand that the data I provided may not seem normal but perhaps it will help if I give some context.

I encountered this issue while developing a sort of “auto ML” platform that automatically selects features and tunes hyperparameters. The feature selection works by starting with one “best” feature, selecting the next best feature, and so on. Hyperparameter tuning is also happening during this process. I ran into the issue during the start of the feature selection process when it was training with only one feature.

If I were doing all of this manually, then yes, maybe it would be reasonable for me to just work around convergence issues. But for an automated platform, I want the underlying learning algorithm implementation (LIBSVM in this case) to be robust enough that it can handle every case, including numerical edge cases.

I hope that makes sense!

[cjlin1]

We have a maximal number of iterations to handle such situations.
Checking such numerical errors first makes the code more complicated
and second slow down the program in the majority of the situations.
Indeed we did not even have a max number of iterations before, and
it was added because of this.
The numerical issues are very rare, so we do not think extra checks are needed.

[fumoboy007]

We have a maximal number of iterations to handle such situations.

Yup, I saw that. However, it’s not great because the training process takes forever to complete. In my case, a regular training process takes less than one second to complete while this problematic one takes hours (I haven’t checked exactly how long).

The numerical issues are very rare, so we do not think extra checks are needed.

Perhaps it’s too early to say whether an extra check is needed? I am trying to investigate the issue more deeply first. Admittedly, it’s difficult for me because I am not an expert here. Would you mind if I ask you some questions about the implementation and my findings?

[cjlin1]

Yes, you are welcome to ask me questions

…

On 2024-12-10 13:25, fumoboy007 wrote: > We have a maximal number of iterations to handle such situations. Yup, I saw that. However, it’s not great because the training process takes forever to complete. In my case, a regular training process takes less than one second to complete while this problematic one takes hours (I haven’t checked exactly how long). > The numerical issues are very rare, so we do not think extra checks > are needed. Perhaps it’s too early to say whether an extra check is needed? I am trying to investigate the issue more deeply first. Admittedly, it’s difficult for me because I am not an expert here. Would you mind if I ask you some questions about the implementation and my findings? -- Reply to this email directly, view it on GitHub [1], or unsubscribe [2]. You are receiving this because you were mentioned.Message ID: ***@***.***> [ { ***@***.***": "http://schema.org", ***@***.***": "EmailMessage", "potentialAction": { ***@***.***": "ViewAction", "target": "#225 (comment)", "url": "#225 (comment)", "name": "View Issue" }, "description": "View this Issue on GitHub", "publisher": { ***@***.***": "Organization", "name": "GitHub", "url": "https://github.com" } } ] Links: ------ [1] #225 (comment) [2] https://github.com/notifications/unsubscribe-auth/ABI3BHU7ZCO6L2K3EPCYU7T2E5L5FAVCNFSM6AAAAABS67EY26VHI2DSMVQWIX3LMV43OSLTON2WKQ3PNVWWK3TUHMZDKMZSHEZDGOBVHE

[fumoboy007]

Thanks! Here are my findings so far. I added some extra logging and focused on the iterations where the training process starts oscillating between i=1525, j=2023 and i=2023, j=1525.

Findings

Iteration 111360

1. Select i=1525 and j=2023 with obj_diff_min=-0.01164691862471009.
2. Calculate quad_coef=0.00000136871642553 and delta=169.13929523467024296.
3. Update alpha and G as follows:

Variable	Previous Value	New Value
alpha[1525]	95.55052799999999991 (upper bound)	34.06510009316370713
alpha[2023]	34.06510009316369292	95.55052799999999991 (upper bound)
G[1525]	1.28270253761288888	3.23039156859958565
G[2023]	1.28247103388129835	3.23109825804289930

Iteration 111361

1. Select i=2023 and j=1525 with obj_diff_min=-0.04086391544522031.
2. Calculate quad_coef=0.00000136871642553 and delta=516.31545448807764842.
3. Update alpha and G as follows:

Variable	Previous Value	New Value
alpha[1525]	34.06510009316370713	95.55052799999999991 (upper bound)
alpha[2023]	95.55052799999999991 (upper bound)	34.06510009316365029
G[1525]	3.23039156859958565	1.28270253760459374
G[2023]	3.23109825804289930	1.28247103387388273

Subsequent Iterations

The values in subsequent iterations are similar to the above two iterations, though the last few digits of obj_diff_min, delta, alpha, and G change slightly.

Questions

• In Iteration 111360, how did both G[1525] and G[2023] increase from 1 to 3? As far as I understand, the gradient for one or both i and j should get closer to the optimal value of 0 after one iteration.
• I mentioned that in subsequent iterations, the last few digits of G change slightly. Those small changes are increases. Again, shouldn’t G get closer to 0, even if very slowly?
• Does the upper bound clipping play any role in the oscillation?
• Why exactly is this oscillation happening? If it is numerical instability, where is it occurring?

[cjlin1]

At optimum, G may not be zero as SVM dual is a constrained
optimization problem. You can see the optimality condition shown
on page 3 of
https://www.csie.ntu.edu.tw/~cjlin/papers/quadworkset.pdf

In theory, things should converge, but apparently numerical
issues occur here. To precisely check what happens, the
best way is to compare with the version that works.
You mentioned that an earlier version stops without issues.
You can print out information in every iteration and do a
comparison.

[fumoboy007]

At optimum, G may not be zero as SVM dual is a constrained optimization problem.

Ohhh I see, that’s why m(α) ≤ M(α) is the stopping condition. Makes sense, thanks!

Will investigate further~~

[fumoboy007]

After a long and painful investigation, I have finally found the root cause of the issue! 🎉

I sent pull request #228 for your consideration. Below is the root cause analysis copied from the pull request description.

Background

The optimization algorithm has three main calculations:

1. Select the working set {i, j} that minimizes the decrease in the objective function.
2. Change alpha[i] and alpha[j] to minimize the decrease in the objective function while respecting constraints.
3. Update the gradient of the objective function according to the changes to alpha[i] and alpha[j].

All three calculations make use of the matrix Q, which is represented by the QMatrix class. The QMatrix class has two main methods:

• get_Q, which returns an array of values for a single column of the matrix; and
• get_QD, which returns an array of diagonal values.

Problem

Q values are of type Qfloat while QD values are of type double. Qfloat is currently defined as float, so there can be inconsistency in the diagonal values returned by get_Q and get_QD. For example, in #225, one of the diagonal values is 181.05748749793070829 as double and 180.99411909539512067 as float.

The first two calculations of the optimization algorithm access the diagonal values via get_QD. However, the third calculation accesses the diagonal values via get_Q. This inconsistency between the minimization calculations and the gradient update can cause the optimization algorithm to oscillate, as demonstrated by #225.

Solution

We change get_Q to return a new class called QColumn instead of a plain array of values. The QColumn class overloads the subscript operator, so accessing individual elements is the same as before. Internally though, the QColumn class will return the QD value when the diagonal element is accessed. This guarantees that all calculations are using the same values for the diagonal elements, eliminating the inconsistency.

Alternatives Considered

Alternatively, we could change Qfloat to be defined as double. This would also eliminate the inconsistency; however, it would reduce the cache capacity by half.