Repository: scikit-learn/scikit-learn
Pull request: [MRG] Expose `n_iter_` to BaseLibSVM
#### Reference Issues/PRs
Fixes #18928

#### What does this implement/fix? Explain your changes.
This PR exposes  `n_iter_` to [BaseLibSVM](https://github.com/scikit-learn/scikit-learn/blob/d4d5f8c7e02cfaced76757fbf38e21c5b28b67b0/sklearn/svm/_base.py#L63):

- Add private attribute `_num_iter` to BaseLibSVM. This will always be a numpy array.
- Add attribute `n_iter_` to BaseLibSVM. When LibSVM optimizes just one model, this is an integer, otherwise, it is a numpy array with the number of iterations for each model.
- Add documentation of `n_iter_` to the classes that inherit from BaseLibSVM (SVR, NuSVR, OneClassSVM, SVC, NuSVC).
- Add `int num_iter` to the class `Solver` in `sklearn/svm/src/libsvm/svm.cpp`. The `Solve` method from this class is responsible for carrying out the optimization. In this method, we save the final number of iterations run in optimization.
- Add `int num_iter` to the struct `decision_function` in `sklearn/svm/src/libsvm/svm.cpp`. Since the svm models do not have direct access to the `Solver` class, we added the attribute `num_iter` to `decision_function`. This struct is returned by the function `svm_train_one` which the train functions call. The function `svm_train_one` calls the Solvers, thus having access to the number of iterations.
- Add `int *num_iter` to `svm_model` and `svm_csr_model` structs in `sklearn/svm/src/libsvm/svm.h`. Since LibSVM can train several models, depending on the model type and the number of classes, here the attribute is an array of integers.
- Fix bug in test `test_libsvm_iris` of `sklearn/svm/tests/test_svm.py`. Since the `fit` function returns a different set of attributes than the ones used by the `predict` function in `sklearn/svm/_libsvm.pyx`, we explicitly selected the correct ones. Note that this didn't fail previously because the `fit_status` variable took the place of the `svm_type` (changing the type of the model).
- Check the attribute `n_iter_` for the classes that inherit from `BaseLibSVM` in `sklearn/utils/estimator_checks.py`. For that I changed the tests `check_non_transformer_estimators_n_iter` and `check_transformer_n_iter`.
- Add changelog
- Test `n_iter_` for all the estimators that inherit from `BaseLibSVM`.

#### Any other comments?

The pre-commit script removed the trailing white spaces in `sklearn/svm/src/libsvm/svm.cpp`.

Problem statement:
Get the number of iterations in SVR
#### Describe the workflow you want to enable
Hi everyone,

I am manipulating `SVR` objects in `GridSearcheCV`. I am able to access the `mean_fit_time` in the `cv_results_`, but I can't access the number of iterations of the optimization problem.

I would like to have this information to properly set the `max_iter` parameter of the `GridSearch`.

#### Describe your proposed solution

I have tried the following:

```python
from sklearn.svm import SVR
from sklearn.datasets import load_boston

# Load data
X, y = load_boston(return_X_y=True)

# Model test
model = SVR(verbose=4)
model.fit(X, y)
```

> [LibSVM]*
optimization finished, #iter = 351
obj = -3012.975812, rho = -21.172739
nSV = 499, nBSV = 431
Out[1]: SVR(gamma=1.0, verbose=4)

I am interested in getting the `#iter` field here. It should be available as a property of the `model` once fitted, and all number of iterations should appear somewhere in the `cv_results_`.

Also, please not that this feature should be available for all libsvm-based SVM objects: `SVC`, `SVR`, etc...

#### Additional context
I am running this code on:
Python 3.7.3
scikit-learn  0.23.1

Thanks by advance for your support.

Hints from issue discussion:
I don't know why we don't expose `n_iter` (it's possible that LibSVM makes it hard for us).

Regarding GridSearch: attributes aren't stored in `cv_results`. We can't store fitted models there because that would mean storing `n_param_candidates` * `n_folds` models, which can be prohibitively expensive in terms of memory.

You can however access the attributes of the final fitted model from `gs.best_estimator_`
For ref:

https://github.com/scikit-learn/scikit-learn/blob/f82525ee2f5ce0d1b014ac5ba83d9d699e299a0b/sklearn/utils/estimator_checks.py#L2751-L2762
> Regarding GridSearch: attributes aren't stored in cv_results. We can't store fitted models there because that would mean storing n_param_candidates * n_folds models, which can be prohibitively expensive in terms of memory.

Actually, we can by passing the parameter `return_estimator=True` and then inspect each estimator.
I am working on an example to show how to do so:

### PR

https://github.com/scikit-learn/scikit-learn/pull/18821

### Artifacts

https://125245-843222-gh.circle-artifacts.com/0/doc/auto_examples/inspection/plot_model_inspection_cross_validation.html
Yup but that's for `cross_validate`, @Leonardbcm mentionned `GridSearchCV` and `cv_resulst_` in particular. `cross_validate` might fit @Leonardbcm 's needs though
Whoops, I misread my bad.

On Fri, 27 Nov 2020 at 20:36, Nicolas Hug <notifications@github.com> wrote:

> Yup but that's for cross_validate, @Leonardbcm
> <https://github.com/Leonardbcm> mentionned GridSearchCV and cv_resulst_
> in particular. cross_validate might fit @Leonardbcm
> <https://github.com/Leonardbcm> 's needs though
>
> —
> You are receiving this because you commented.
> Reply to this email directly, view it on GitHub
> <https://github.com/scikit-learn/scikit-learn/issues/18928#issuecomment-734962975>,
> or unsubscribe
> <https://github.com/notifications/unsubscribe-auth/ABY32PZ6XQ7YK3QWF2GV34LSR75URANCNFSM4UEUNKWA>
> .
>


-- 
Guillaume Lemaitre
Scikit-learn @ Inria Foundation
https://glemaitre.github.io/

Hacky solution: add a scorer under the grid searcher's scoring parameter, defined like:
```py
def get_n_iter(est, X_test, y_test):
    return est.n_iter_
```
Hi, thanks for your answers.

> I don't know why we don't expose n_iter (it's possible that LibSVM makes it hard for us).

This is sad but it makes a lot of sense.

> Regarding GridSearch: attributes aren't stored in cv_results. We can't store fitted models there because that would mean storing n_param_candidates * n_folds models, which can be prohibitively expensive in terms of memory.

It also makes sense but could be optional. Metrics values are already stored and they also take `n_metrics` * `n_folds` per grid search iteration.

> cross_validate might fit @Leonardbcm 's needs though

Unfortunately, no. My use-case is the following: I run my grid search on 20CPUs. Most of the iterations are taking 10s per fold to train and evaluate, while a few other are taking hours to converge. As a single model training can't be dispatched on multiple cores, the result is that I have 18 CPUs that are done and 2 that "blocks" the main script. (Does that make sense?).

What I am trying to do is setting the `max_iter` hyper-parameter for my entiere grid search so that the search is not "blocked" by models that are hard to converge. Accessing the number of realized itertions would allow me to:
- FInd a reasonnable value for `max_iter` to allow most model to converge while preserving a low training time
- Know which models have finished training but have not reached convergence yet (their number of iteration is equal to `max_iter`.

 > return est.n_iter_

The problem is that `n_iter_` is not an attribute nor a property of the models.

Thanks for your help.



If it is ok, I would like to work on this.
I already have a draft with the number of iterations stored as an attribute for the class `BaseLibSVM`.
Though, I still need to document the changes and do more testing.

Should we try to expose more attributes from the libsvm code-base?

Exposing `n_iter_` would be enough for the moment.

Commit message: Comment unpacking and repacking of values for test

Co-authored-by: Julien Jerphanion <git@jjerphan.xyz>

Patch to review:
diff --git a/doc/whats_new/v1.1.rst b/doc/whats_new/v1.1.rst
index 2cf9d9226dbff..66e0e7c6610a0 100644
--- a/doc/whats_new/v1.1.rst
+++ b/doc/whats_new/v1.1.rst
@@ -124,6 +124,14 @@ Changelog
   Setting a transformer to "passthrough" will pass the features unchanged.
   :pr:`20860` by :user:`Shubhraneel Pal <shubhraneel>`.
 
+:mod:`sklearn.svm`
+...................
+
+- |Enhancement| :class:`svm.OneClassSVM`, :class:`svm.NuSVC`,
+  :class:`svm.NuSVR`, :class:`svm.SVC` and :class:`svm.SVR` now expose
+  `n_iter_`, the number of iterations used by the libsvm optimization routine.
+  :pr:`21408` by :user:`Juan Martín Loyola <jmloyola>`.
+
 :mod:`sklearn.utils`
 ....................
 
diff --git a/sklearn/svm/_base.py b/sklearn/svm/_base.py
index 7f0e628aab70c..614447e4467a1 100644
--- a/sklearn/svm/_base.py
+++ b/sklearn/svm/_base.py
@@ -266,6 +266,15 @@ def fit(self, X, y, sample_weight=None):
             self.intercept_ *= -1
             self.dual_coef_ = -self.dual_coef_
 
+        # If the number of models optimized by libSVM is one, get the number of
+        # iterations as an integer instead of ndarray.
+        if (
+            self._impl in ["c_svc", "nu_svc"] and len(self.classes_) <= 2
+        ) or self._impl in ["one_class", "epsilon_svr", "nu_svr"]:
+            self.n_iter_ = self._num_iter[0]
+        else:
+            self.n_iter_ = self._num_iter
+
         return self
 
     def _validate_targets(self, y):
@@ -312,6 +321,7 @@ def _dense_fit(self, X, y, sample_weight, solver_type, kernel, random_seed):
             self._probA,
             self._probB,
             self.fit_status_,
+            self._num_iter,
         ) = libsvm.fit(
             X,
             y,
@@ -352,6 +362,7 @@ def _sparse_fit(self, X, y, sample_weight, solver_type, kernel, random_seed):
             self._probA,
             self._probB,
             self.fit_status_,
+            self._num_iter,
         ) = libsvm_sparse.libsvm_sparse_train(
             X.shape[1],
             X.data,
@@ -439,6 +450,7 @@ def _dense_predict(self, X):
             self._n_support,
             self._dual_coef_,
             self._intercept_,
+            self._num_iter,
             self._probA,
             self._probB,
             svm_type=svm_type,
@@ -483,6 +495,7 @@ def _sparse_predict(self, X):
             self._n_support,
             self._probA,
             self._probB,
+            self._num_iter,
         )
 
     def _compute_kernel(self, X):
@@ -540,6 +553,7 @@ def _dense_decision_function(self, X):
             self._n_support,
             self._dual_coef_,
             self._intercept_,
+            self._num_iter,
             self._probA,
             self._probB,
             svm_type=LIBSVM_IMPL.index(self._impl),
@@ -583,6 +597,7 @@ def _sparse_decision_function(self, X):
             self._n_support,
             self._probA,
             self._probB,
+            self._num_iter,
         )
 
     def _validate_for_predict(self, X):
@@ -886,6 +901,7 @@ def _dense_predict_proba(self, X):
             self._n_support,
             self._dual_coef_,
             self._intercept_,
+            self._num_iter,
             self._probA,
             self._probB,
             svm_type=svm_type,
@@ -931,6 +947,7 @@ def _sparse_predict_proba(self, X):
             self._n_support,
             self._probA,
             self._probB,
+            self._num_iter,
         )
 
     def _get_coef(self):
diff --git a/sklearn/svm/_classes.py b/sklearn/svm/_classes.py
index 719540cd725c0..4f4bce25f88c3 100644
--- a/sklearn/svm/_classes.py
+++ b/sklearn/svm/_classes.py
@@ -670,6 +670,15 @@ class SVC(BaseSVC):
 
         .. versionadded:: 1.0
 
+    n_iter_ : int or ndarray of shape (n_class * (n_class - 1) // 2,)
+        Number of iterations run by the optimization routine to fit the model.
+        If `classes_ <= 2`, only one model is optimized, thus an integer is
+        returned. Otherwise, multiple models are optimized separately, thus
+        having multiple number of iterations. In this case a numpy array is
+        returned with the number of iterations for each model.
+
+        .. versionadded:: 1.1
+
     support_ : ndarray of shape (n_SV)
         Indices of support vectors.
 
@@ -925,6 +934,15 @@ class NuSVC(BaseSVC):
 
         .. versionadded:: 1.0
 
+    n_iter_ : int or ndarray of shape (n_class * (n_class - 1) // 2,)
+        Number of iterations run by the optimization routine to fit the model.
+        If `classes_ <= 2`, only one model is optimized, thus an integer is
+        returned. Otherwise, multiple models are optimized separately, thus
+        having multiple number of iterations. In this case a numpy array is
+        returned with the number of iterations for each model.
+
+        .. versionadded:: 1.1
+
     support_ : ndarray of shape (n_SV,)
         Indices of support vectors.
 
@@ -1140,6 +1158,11 @@ class SVR(RegressorMixin, BaseLibSVM):
 
         .. versionadded:: 1.0
 
+    n_iter_ : int
+        Number of iterations run by the optimization routine to fit the model.
+
+        .. versionadded:: 1.1
+
     n_support_ : ndarray of shape (n_classes,), dtype=int32
         Number of support vectors for each class.
 
@@ -1328,6 +1351,11 @@ class NuSVR(RegressorMixin, BaseLibSVM):
 
         .. versionadded:: 1.0
 
+    n_iter_ : int
+        Number of iterations run by the optimization routine to fit the model.
+
+        .. versionadded:: 1.1
+
     n_support_ : ndarray of shape (n_classes,), dtype=int32
         Number of support vectors for each class.
 
@@ -1512,6 +1540,11 @@ class OneClassSVM(OutlierMixin, BaseLibSVM):
 
         .. versionadded:: 1.0
 
+    n_iter_ : int
+        Number of iterations run by the optimization routine to fit the model.
+
+        .. versionadded:: 1.1
+
     n_support_ : ndarray of shape (n_classes,), dtype=int32
         Number of support vectors for each class.
 
diff --git a/sklearn/svm/_libsvm.pxi b/sklearn/svm/_libsvm.pxi
index ab7d212e3ba1a..39cc40a0bc9ae 100644
--- a/sklearn/svm/_libsvm.pxi
+++ b/sklearn/svm/_libsvm.pxi
@@ -53,9 +53,10 @@ cdef extern from "libsvm_helper.c":
 
     svm_model *set_model (svm_parameter *, int, char *, np.npy_intp *,
                          char *, np.npy_intp *, np.npy_intp *, char *,
-                         char *, char *, char *, char *)
+                         char *, char *, char *, char *, char *)
 
     void copy_sv_coef   (char *, svm_model *)
+    void copy_n_iter  (char *, svm_model *)
     void copy_intercept (char *, svm_model *, np.npy_intp *)
     void copy_SV        (char *, svm_model *, np.npy_intp *)
     int copy_support (char *data, svm_model *model)
diff --git a/sklearn/svm/_libsvm.pyx b/sklearn/svm/_libsvm.pyx
index 9488bda4ccf58..cbe1fe9893564 100644
--- a/sklearn/svm/_libsvm.pyx
+++ b/sklearn/svm/_libsvm.pyx
@@ -155,6 +155,9 @@ def fit(
 
     probA, probB : array of shape (n_class*(n_class-1)/2,)
         Probability estimates, empty array for probability=False.
+
+    n_iter : ndarray of shape (max(1, (n_class * (n_class - 1) // 2)),)
+        Number of iterations run by the optimization routine to fit the model.
     """
 
     cdef svm_parameter param
@@ -202,6 +205,10 @@ def fit(
     SV_len  = get_l(model)
     n_class = get_nr(model)
 
+    cdef np.ndarray[int, ndim=1, mode='c'] n_iter
+    n_iter = np.empty(max(1, n_class * (n_class - 1) // 2), dtype=np.intc)
+    copy_n_iter(n_iter.data, model)
+
     cdef np.ndarray[np.float64_t, ndim=2, mode='c'] sv_coef
     sv_coef = np.empty((n_class-1, SV_len), dtype=np.float64)
     copy_sv_coef (sv_coef.data, model)
@@ -251,7 +258,7 @@ def fit(
     free(problem.x)
 
     return (support, support_vectors, n_class_SV, sv_coef, intercept,
-           probA, probB, fit_status)
+           probA, probB, fit_status, n_iter)
 
 
 cdef void set_predict_params(
@@ -282,6 +289,7 @@ def predict(np.ndarray[np.float64_t, ndim=2, mode='c'] X,
             np.ndarray[np.int32_t, ndim=1, mode='c'] nSV,
             np.ndarray[np.float64_t, ndim=2, mode='c'] sv_coef,
             np.ndarray[np.float64_t, ndim=1, mode='c'] intercept,
+            np.ndarray[int, ndim=1, mode='c'] n_iter,
             np.ndarray[np.float64_t, ndim=1, mode='c'] probA=np.empty(0),
             np.ndarray[np.float64_t, ndim=1, mode='c'] probB=np.empty(0),
             int svm_type=0, kernel='rbf', int degree=3,
@@ -313,6 +321,9 @@ def predict(np.ndarray[np.float64_t, ndim=2, mode='c'] X,
     intercept : array of shape (n_class*(n_class-1)/2)
         Intercept in decision function.
 
+    n_iter : ndarray of shape (max(1, (n_class * (n_class - 1) // 2)),)
+        Number of iterations run by the optimization routine to fit the model.
+
     probA, probB : array of shape (n_class*(n_class-1)/2,)
         Probability estimates.
 
@@ -353,7 +364,8 @@ def predict(np.ndarray[np.float64_t, ndim=2, mode='c'] X,
                        class_weight_label.data, class_weight.data)
     model = set_model(&param, <int> nSV.shape[0], SV.data, SV.shape,
                       support.data, support.shape, sv_coef.strides,
-                      sv_coef.data, intercept.data, nSV.data, probA.data, probB.data)
+                      sv_coef.data, intercept.data, nSV.data, probA.data,
+                      probB.data, n_iter.data)
     cdef BlasFunctions blas_functions
     blas_functions.dot = _dot[double]
     #TODO: use check_model
@@ -376,6 +388,7 @@ def predict_proba(
     np.ndarray[np.int32_t, ndim=1, mode='c'] nSV,
     np.ndarray[np.float64_t, ndim=2, mode='c'] sv_coef,
     np.ndarray[np.float64_t, ndim=1, mode='c'] intercept,
+    np.ndarray[int, ndim=1, mode='c'] n_iter,
     np.ndarray[np.float64_t, ndim=1, mode='c'] probA=np.empty(0),
     np.ndarray[np.float64_t, ndim=1, mode='c'] probB=np.empty(0),
     int svm_type=0, kernel='rbf', int degree=3,
@@ -417,6 +430,9 @@ def predict_proba(
     intercept : array of shape (n_class*(n_class-1)/2,)
         Intercept in decision function.
 
+    n_iter : ndarray of shape (max(1, (n_class * (n_class - 1) // 2)),)
+        Number of iterations run by the optimization routine to fit the model.
+
     probA, probB : array of shape (n_class*(n_class-1)/2,)
         Probability estimates.
 
@@ -457,7 +473,7 @@ def predict_proba(
     model = set_model(&param, <int> nSV.shape[0], SV.data, SV.shape,
                       support.data, support.shape, sv_coef.strides,
                       sv_coef.data, intercept.data, nSV.data,
-                      probA.data, probB.data)
+                      probA.data, probB.data, n_iter.data)
 
     cdef np.npy_intp n_class = get_nr(model)
     cdef BlasFunctions blas_functions
@@ -481,6 +497,7 @@ def decision_function(
     np.ndarray[np.int32_t, ndim=1, mode='c'] nSV,
     np.ndarray[np.float64_t, ndim=2, mode='c'] sv_coef,
     np.ndarray[np.float64_t, ndim=1, mode='c'] intercept,
+    np.ndarray[int, ndim=1, mode='c'] n_iter,
     np.ndarray[np.float64_t, ndim=1, mode='c'] probA=np.empty(0),
     np.ndarray[np.float64_t, ndim=1, mode='c'] probB=np.empty(0),
     int svm_type=0, kernel='rbf', int degree=3,
@@ -515,6 +532,9 @@ def decision_function(
     intercept : array, shape=[n_class*(n_class-1)/2]
         Intercept in decision function.
 
+    n_iter : ndarray of shape (max(1, (n_class * (n_class - 1) // 2)),)
+        Number of iterations run by the optimization routine to fit the model.
+
     probA, probB : array, shape=[n_class*(n_class-1)/2]
         Probability estimates.
 
@@ -559,7 +579,7 @@ def decision_function(
     model = set_model(&param, <int> nSV.shape[0], SV.data, SV.shape,
                       support.data, support.shape, sv_coef.strides,
                       sv_coef.data, intercept.data, nSV.data,
-                      probA.data, probB.data)
+                      probA.data, probB.data, n_iter.data)
 
     if svm_type > 1:
         n_class = 1
diff --git a/sklearn/svm/_libsvm_sparse.pyx b/sklearn/svm/_libsvm_sparse.pyx
index 92d94b0c685a5..b73c2cd487a4b 100644
--- a/sklearn/svm/_libsvm_sparse.pyx
+++ b/sklearn/svm/_libsvm_sparse.pyx
@@ -35,12 +35,13 @@ cdef extern from "libsvm_sparse_helper.c":
                             char *SV_indices, np.npy_intp *SV_intptr_dims,
                             char *SV_intptr,
                             char *sv_coef, char *rho, char *nSV,
-                            char *probA, char *probB)
+                            char *probA, char *probB, char *n_iter)
     svm_parameter *set_parameter (int , int , int , double, double ,
                                   double , double , double , double,
                                   double, int, int, int, char *, char *, int,
                                   int)
     void copy_sv_coef   (char *, svm_csr_model *)
+    void copy_n_iter  (char *, svm_csr_model *)
     void copy_support   (char *, svm_csr_model *)
     void copy_intercept (char *, svm_csr_model *, np.npy_intp *)
     int copy_predict (char *, svm_csr_model *, np.npy_intp *, char *, BlasFunctions *)
@@ -159,6 +160,10 @@ def libsvm_sparse_train ( int n_features,
     cdef np.npy_intp SV_len = get_l(model)
     cdef np.npy_intp n_class = get_nr(model)
 
+    cdef np.ndarray[int, ndim=1, mode='c'] n_iter
+    n_iter = np.empty(max(1, n_class * (n_class - 1) // 2), dtype=np.intc)
+    copy_n_iter(n_iter.data, model)
+
     # copy model.sv_coef
     # we create a new array instead of resizing, otherwise
     # it would not erase previous information
@@ -217,7 +222,7 @@ def libsvm_sparse_train ( int n_features,
     free_param(param)
 
     return (support, support_vectors_, sv_coef_data, intercept, n_class_SV,
-            probA, probB, fit_status)
+            probA, probB, fit_status, n_iter)
 
 
 def libsvm_sparse_predict (np.ndarray[np.float64_t, ndim=1, mode='c'] T_data,
@@ -236,7 +241,8 @@ def libsvm_sparse_predict (np.ndarray[np.float64_t, ndim=1, mode='c'] T_data,
                             shrinking, int probability,
                             np.ndarray[np.int32_t, ndim=1, mode='c'] nSV,
                             np.ndarray[np.float64_t, ndim=1, mode='c'] probA,
-                            np.ndarray[np.float64_t, ndim=1, mode='c'] probB):
+                            np.ndarray[np.float64_t, ndim=1, mode='c'] probB,
+                            np.ndarray[int, ndim=1, mode='c'] n_iter):
     """
     Predict values T given a model.
 
@@ -277,7 +283,7 @@ def libsvm_sparse_predict (np.ndarray[np.float64_t, ndim=1, mode='c'] T_data,
                           SV_indices.shape, SV_indices.data,
                           SV_indptr.shape, SV_indptr.data,
                           sv_coef.data, intercept.data,
-                          nSV.data, probA.data, probB.data)
+                          nSV.data, probA.data, probB.data, n_iter.data)
     #TODO: use check_model
     dec_values = np.empty(T_indptr.shape[0]-1)
     cdef BlasFunctions blas_functions
@@ -313,7 +319,8 @@ def libsvm_sparse_predict_proba(
     double nu, double p, int shrinking, int probability,
     np.ndarray[np.int32_t, ndim=1, mode='c'] nSV,
     np.ndarray[np.float64_t, ndim=1, mode='c'] probA,
-    np.ndarray[np.float64_t, ndim=1, mode='c'] probB):
+    np.ndarray[np.float64_t, ndim=1, mode='c'] probB,
+    np.ndarray[int, ndim=1, mode='c'] n_iter):
     """
     Predict values T given a model.
     """
@@ -334,7 +341,7 @@ def libsvm_sparse_predict_proba(
                           SV_indices.shape, SV_indices.data,
                           SV_indptr.shape, SV_indptr.data,
                           sv_coef.data, intercept.data,
-                          nSV.data, probA.data, probB.data)
+                          nSV.data, probA.data, probB.data, n_iter.data)
     #TODO: use check_model
     cdef np.npy_intp n_class = get_nr(model)
     cdef int rv
@@ -374,7 +381,8 @@ def libsvm_sparse_decision_function(
     double nu, double p, int shrinking, int probability,
     np.ndarray[np.int32_t, ndim=1, mode='c'] nSV,
     np.ndarray[np.float64_t, ndim=1, mode='c'] probA,
-    np.ndarray[np.float64_t, ndim=1, mode='c'] probB):
+    np.ndarray[np.float64_t, ndim=1, mode='c'] probB,
+    np.ndarray[int, ndim=1, mode='c'] n_iter):
     """
     Predict margin (libsvm name for this is predict_values)
 
@@ -399,7 +407,7 @@ def libsvm_sparse_decision_function(
                           SV_indices.shape, SV_indices.data,
                           SV_indptr.shape, SV_indptr.data,
                           sv_coef.data, intercept.data,
-                          nSV.data, probA.data, probB.data)
+                          nSV.data, probA.data, probB.data, n_iter.data)
 
     if svm_type > 1:
         n_class = 1
diff --git a/sklearn/svm/src/libsvm/LIBSVM_CHANGES b/sklearn/svm/src/libsvm/LIBSVM_CHANGES
index bde6beaca2694..663550b8ddd6f 100644
--- a/sklearn/svm/src/libsvm/LIBSVM_CHANGES
+++ b/sklearn/svm/src/libsvm/LIBSVM_CHANGES
@@ -7,4 +7,5 @@ This is here mainly as checklist for incorporation of new versions of libsvm.
   * Improved random number generator (fix on windows, enhancement on other
     platforms). See <https://github.com/scikit-learn/scikit-learn/pull/13511#issuecomment-481729756>
   * invoke scipy blas api for svm kernel function to improve performance with speedup rate of 1.5X to 2X for dense data only. See <https://github.com/scikit-learn/scikit-learn/pull/16530>
+  * Expose the number of iterations run in optimization. See <https://github.com/scikit-learn/scikit-learn/pull/21408>
 The changes made with respect to upstream are detailed in the heading of svm.cpp
diff --git a/sklearn/svm/src/libsvm/libsvm_helper.c b/sklearn/svm/src/libsvm/libsvm_helper.c
index 17f328f9e7c4c..c11ff6703fd77 100644
--- a/sklearn/svm/src/libsvm/libsvm_helper.c
+++ b/sklearn/svm/src/libsvm/libsvm_helper.c
@@ -2,6 +2,13 @@
 #include <numpy/arrayobject.h>
 #include "svm.h"
 #include "_svm_cython_blas_helpers.h"
+
+
+#ifndef max
+    #define max(x, y) (((x) > (y)) ? (x) : (y))
+#endif
+
+
 /*
  * Some helper methods for libsvm bindings.
  *
@@ -109,13 +116,14 @@ struct svm_model *set_model(struct svm_parameter *param, int nr_class,
                             char *support, npy_intp *support_dims,
                             npy_intp *sv_coef_strides,
                             char *sv_coef, char *rho, char *nSV,
-                            char *probA, char *probB)
+                            char *probA, char *probB, char *n_iter)
 {
     struct svm_model *model;
     double *dsv_coef = (double *) sv_coef;
-    int i, m;
+    int i, m, n_models;
 
     m = nr_class * (nr_class-1)/2;
+    n_models = nr_class <= 2 ? 1 : m;
 
     if ((model = malloc(sizeof(struct svm_model))) == NULL)
         goto model_error;
@@ -127,6 +135,8 @@ struct svm_model *set_model(struct svm_parameter *param, int nr_class,
         goto sv_coef_error;
     if ((model->rho = malloc( m * sizeof(double))) == NULL)
         goto rho_error;
+    if ((model->n_iter = malloc(n_models * sizeof(int))) == NULL)
+        goto n_iter_error;
 
     model->nr_class = nr_class;
     model->param = *param;
@@ -181,6 +191,8 @@ struct svm_model *set_model(struct svm_parameter *param, int nr_class,
     model->free_sv = 0;
     return model;
 
+n_iter_error:
+    free(model->n_iter);
 probB_error:
     free(model->probA);
 probA_error:
@@ -218,6 +230,15 @@ npy_intp get_nr(struct svm_model *model)
     return (npy_intp) model->nr_class;
 }
 
+/*
+ * Get the number of iterations run in optimization
+ */
+void copy_n_iter(char *data, struct svm_model *model)
+{
+    const int n_models = max(1, model->nr_class * (model->nr_class-1)/2);
+    memcpy(data, model->n_iter, n_models * sizeof(int));
+}
+
 /*
  * Some helpers to convert from libsvm sparse data structures
  * model->sv_coef is a double **, whereas data is just a double *,
@@ -371,6 +392,7 @@ int free_model(struct svm_model *model)
     free(model->label);
     free(model->probA);
     free(model->probB);
+    free(model->n_iter);
     free(model->nSV);
     free(model);
 
diff --git a/sklearn/svm/src/libsvm/libsvm_sparse_helper.c b/sklearn/svm/src/libsvm/libsvm_sparse_helper.c
index a85a532319d88..e94dc350ff3c2 100644
--- a/sklearn/svm/src/libsvm/libsvm_sparse_helper.c
+++ b/sklearn/svm/src/libsvm/libsvm_sparse_helper.c
@@ -3,6 +3,12 @@
 #include "svm.h"
 #include "_svm_cython_blas_helpers.h"
 
+
+#ifndef max
+    #define max(x, y) (((x) > (y)) ? (x) : (y))
+#endif
+
+
 /*
  * Convert scipy.sparse.csr to libsvm's sparse data structure
  */
@@ -103,13 +109,14 @@ struct svm_csr_model *csr_set_model(struct svm_parameter *param, int nr_class,
                             char *SV_indices, npy_intp *SV_indptr_dims,
                             char *SV_intptr,
                             char *sv_coef, char *rho, char *nSV,
-                            char *probA, char *probB)
+                            char *probA, char *probB, char *n_iter)
 {
     struct svm_csr_model *model;
     double *dsv_coef = (double *) sv_coef;
-    int i, m;
+    int i, m, n_models;
 
     m = nr_class * (nr_class-1)/2;
+    n_models = nr_class <= 2 ? 1 : m;
 
     if ((model = malloc(sizeof(struct svm_csr_model))) == NULL)
         goto model_error;
@@ -121,6 +128,8 @@ struct svm_csr_model *csr_set_model(struct svm_parameter *param, int nr_class,
         goto sv_coef_error;
     if ((model->rho = malloc( m * sizeof(double))) == NULL)
         goto rho_error;
+    if ((model->n_iter = malloc(n_models * sizeof(int))) == NULL)
+        goto n_iter_error;
 
     /* in the case of precomputed kernels we do not use
        dense_to_precomputed because we don't want the leading 0. As
@@ -180,6 +189,8 @@ struct svm_csr_model *csr_set_model(struct svm_parameter *param, int nr_class,
     model->free_sv = 0;
     return model;
 
+n_iter_error:
+    free(model->n_iter);
 probB_error:
     free(model->probA);
 probA_error:
@@ -348,6 +359,15 @@ void copy_sv_coef(char *data, struct svm_csr_model *model)
     }
 }
 
+/*
+ * Get the number of iterations run in optimization
+ */
+void copy_n_iter(char *data, struct svm_csr_model *model)
+{
+    const int n_models = max(1, model->nr_class * (model->nr_class-1)/2);
+    memcpy(data, model->n_iter, n_models * sizeof(int));
+}
+
 /*
  * Get the number of support vectors in a model.
  */
@@ -409,6 +429,7 @@ int free_model(struct svm_csr_model *model)
     free(model->label);
     free(model->probA);
     free(model->probB);
+    free(model->n_iter);
     free(model->nSV);
     free(model);
 
diff --git a/sklearn/svm/src/libsvm/svm.cpp b/sklearn/svm/src/libsvm/svm.cpp
index d209e35fc0a35..de07fecdba2ac 100644
--- a/sklearn/svm/src/libsvm/svm.cpp
+++ b/sklearn/svm/src/libsvm/svm.cpp
@@ -55,6 +55,10 @@ SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
      Sylvain Marie, Schneider Electric
      see <https://github.com/scikit-learn/scikit-learn/pull/13511#issuecomment-481729756>
 
+   Modified 2021:
+
+   - Exposed number of iterations run in optimization, Juan Martín Loyola.
+     See <https://github.com/scikit-learn/scikit-learn/pull/21408/>
  */
 
 #include <math.h>
@@ -553,6 +557,7 @@ class Solver {
                 double *upper_bound;
 		double r;	// for Solver_NU
                 bool solve_timed_out;
+		int n_iter;
 	};
 
 	void Solve(int l, const QMatrix& Q, const double *p_, const schar *y_,
@@ -919,6 +924,9 @@ void Solver::Solve(int l, const QMatrix& Q, const double *p_, const schar *y_,
 	for(int i=0;i<l;i++)
 		si->upper_bound[i] = C[i];
 
+	// store number of iterations
+	si->n_iter = iter;
+
 	info("\noptimization finished, #iter = %d\n",iter);
 
 	delete[] p;
@@ -1837,6 +1845,7 @@ struct decision_function
 {
 	double *alpha;
 	double rho;	
+	int n_iter;
 };
 
 static decision_function svm_train_one(
@@ -1902,6 +1911,7 @@ static decision_function svm_train_one(
 	decision_function f;
 	f.alpha = alpha;
 	f.rho = si.rho;
+	f.n_iter = si.n_iter;
 	return f;
 }
 
@@ -2387,6 +2397,8 @@ PREFIX(model) *PREFIX(train)(const PREFIX(problem) *prob, const svm_parameter *p
                 NAMESPACE::decision_function f = NAMESPACE::svm_train_one(prob,param,0,0, status,blas_functions);
 		model->rho = Malloc(double,1);
 		model->rho[0] = f.rho;
+		model->n_iter = Malloc(int,1);
+		model->n_iter[0] = f.n_iter;
 
 		int nSV = 0;
 		int i;
@@ -2523,8 +2535,12 @@ PREFIX(model) *PREFIX(train)(const PREFIX(problem) *prob, const svm_parameter *p
 			model->label[i] = label[i];
 		
 		model->rho = Malloc(double,nr_class*(nr_class-1)/2);
+		model->n_iter = Malloc(int,nr_class*(nr_class-1)/2);
 		for(i=0;i<nr_class*(nr_class-1)/2;i++)
+		{
 			model->rho[i] = f[i].rho;
+			model->n_iter[i] = f[i].n_iter;
+		}
 
 		if(param->probability)
 		{
@@ -2978,6 +2994,9 @@ void PREFIX(free_model_content)(PREFIX(model)* model_ptr)
 
 	free(model_ptr->nSV);
 	model_ptr->nSV = NULL;
+
+	free(model_ptr->n_iter);
+	model_ptr->n_iter = NULL;
 }
 
 void PREFIX(free_and_destroy_model)(PREFIX(model)** model_ptr_ptr)
diff --git a/sklearn/svm/src/libsvm/svm.h b/sklearn/svm/src/libsvm/svm.h
index a1634119858f1..518872c67bc5c 100644
--- a/sklearn/svm/src/libsvm/svm.h
+++ b/sklearn/svm/src/libsvm/svm.h
@@ -76,6 +76,7 @@ struct svm_model
 	int l;			/* total #SV */
 	struct svm_node *SV;		/* SVs (SV[l]) */
 	double **sv_coef;	/* coefficients for SVs in decision functions (sv_coef[k-1][l]) */
+	int *n_iter;		/* number of iterations run by the optimization routine to fit the model */
 
 	int *sv_ind;            /* index of support vectors */
 
@@ -101,6 +102,7 @@ struct svm_csr_model
 	int l;			/* total #SV */
 	struct svm_csr_node **SV;		/* SVs (SV[l]) */
 	double **sv_coef;	/* coefficients for SVs in decision functions (sv_coef[k-1][l]) */
+	int *n_iter;		/* number of iterations run by the optimization routine to fit the model */
 
         int *sv_ind;            /* index of support vectors */
 
diff --git a/sklearn/svm/tests/test_svm.py b/sklearn/svm/tests/test_svm.py
index eb15aca0096b8..a6b77a4ef98d6 100644
--- a/sklearn/svm/tests/test_svm.py
+++ b/sklearn/svm/tests/test_svm.py
@@ -66,12 +66,60 @@ def test_libsvm_iris():
     assert_array_equal(clf.classes_, np.sort(clf.classes_))
 
     # check also the low-level API
-    model = _libsvm.fit(iris.data, iris.target.astype(np.float64))
-    pred = _libsvm.predict(iris.data, *model)
+    # We unpack the values to be able to repack some of the
+    # return values from libsvm.
+    (
+        libsvm_support,
+        libsvm_support_vectors,
+        libsvm_n_class_SV,
+        libsvm_sv_coef,
+        libsvm_intercept,
+        libsvm_probA,
+        libsvm_probB,
+        # libsvm_fit_status won't be packed bellow.
+        libsvm_fit_status,
+        libsvm_n_iter,
+    ) = _libsvm.fit(iris.data, iris.target.astype(np.float64))
+
+    lib_svm_model = (
+        libsvm_support,
+        libsvm_support_vectors,
+        libsvm_n_class_SV,
+        libsvm_sv_coef,
+        libsvm_intercept,
+        libsvm_n_iter,
+        libsvm_probA,
+        libsvm_probB,
+    )
+    pred = _libsvm.predict(iris.data, *lib_svm_model)
     assert np.mean(pred == iris.target) > 0.95
 
-    model = _libsvm.fit(iris.data, iris.target.astype(np.float64), kernel="linear")
-    pred = _libsvm.predict(iris.data, *model, kernel="linear")
+    # We unpack the values to be able to repack some of the
+    # return values from libsvm.
+    (
+        libsvm_support,
+        libsvm_support_vectors,
+        libsvm_n_class_SV,
+        libsvm_sv_coef,
+        libsvm_intercept,
+        libsvm_probA,
+        libsvm_probB,
+        # libsvm_fit_status won't be packed bellow.
+        libsvm_fit_status,
+        libsvm_n_iter,
+    ) = _libsvm.fit(iris.data, iris.target.astype(np.float64), kernel="linear")
+
+    lib_svm_model = (
+        libsvm_support,
+        libsvm_support_vectors,
+        libsvm_n_class_SV,
+        libsvm_sv_coef,
+        libsvm_intercept,
+        libsvm_n_iter,
+        libsvm_probA,
+        libsvm_probB,
+    )
+    pred = _libsvm.predict(iris.data, *lib_svm_model, kernel="linear")
     assert np.mean(pred == iris.target) > 0.95
 
     pred = _libsvm.cross_validation(
diff --git a/sklearn/utils/estimator_checks.py b/sklearn/utils/estimator_checks.py
index 4eb9868d4c9cf..94d8010ec1e82 100644
--- a/sklearn/utils/estimator_checks.py
+++ b/sklearn/utils/estimator_checks.py
@@ -3239,11 +3239,7 @@ def check_non_transformer_estimators_n_iter(name, estimator_orig):
     # labeled, hence n_iter_ = 0 is valid.
     not_run_check_n_iter = [
         "Ridge",
-        "SVR",
-        "NuSVR",
-        "NuSVC",
         "RidgeClassifier",
-        "SVC",
         "RandomizedLasso",
         "LogisticRegressionCV",
         "LinearSVC",
@@ -3268,9 +3264,16 @@ def check_non_transformer_estimators_n_iter(name, estimator_orig):
 
         set_random_state(estimator, 0)
 
+        if name in ["SVC", "NuSVC"] and estimator.kernel == "precomputed":
+            X = np.dot(X, X.T)
+
         estimator.fit(X, y_)
 
-        assert estimator.n_iter_ >= 1
+        # These return a n_iter per model optimized
+        if name in ["SVC", "NuSVC"] and len(estimator.classes_) > 2:
+            assert np.all(estimator.n_iter_ >= 1)
+        else:
+            assert estimator.n_iter_ >= 1
 
 
 @ignore_warnings(category=FutureWarning)
@@ -3297,7 +3300,9 @@ def check_transformer_n_iter(name, estimator_orig):
         estimator.fit(X, y_)
 
         # These return a n_iter per component.
-        if name in CROSS_DECOMPOSITION:
+        if (name in CROSS_DECOMPOSITION) or (
+            name in ["SVC", "NuSVC"] and len(estimator.classes_) > 2
+        ):
             for iter_ in estimator.n_iter_:
                 assert iter_ >= 1
         else:


Code review comments:
diff_hunk: @@ -3268,9 +3264,16 @@ def check_non_transformer_estimators_n_iter(name, estimator_orig):
 
         set_random_state(estimator, 0)
 
+        if name in ["SVC", "NuSVC"] and estimator.kernel == "precomputed":
+            X = np.dot(X, X.T)
path: sklearn/utils/estimator_checks.py
text: @user1:
We have a dedicated `_pairwise_estimator_convert_X` utility for this:

```suggestion
        X = _pairwise_estimator_convert_X(X, estimator_orig)
```

This way, there is no need to hard-code the name of the estimator in the check and it will also work for `SVR` / `NuSVR`.

diff_hunk: @@ -3,6 +3,12 @@
 #include "svm.h"
 #include "_svm_cython_blas_helpers.h"
 
+
+#ifndef max
+    #define max(x, y) (((x) > (y)) ? (x) : (y))
path: sklearn/svm/src/libsvm/libsvm_sparse_helper.c
text: @user1:
Are the inner parenthesis really necessary? I find the following easier to parse:

```suggestion
    #define max(x, y) ((x > y) ? x : y)
```

@user1:
see also: https://github.com/scikit-learn/scikit-learn/pull/21408/files#r740018186

@user2:
They are needed for correct operator precedence.

For instance, it is better that
```c
max(1+1, 2+2)
```
gets expended to:
```c
(((1+1) > (2+2)) ? (1+1) : (2+2))
```

Though this is a simple example and this might not cause problem for this macro.

@author:
As @user2 mentioned, they are needed for correct operator precedence.
In case one of `x` or `y` has an operator with less precedence than the relational operators (<, >), the result will not be the expected one.
Though, it is true that in this case, this does not happen and probably will never happen (unless someone calls the macro with something like `max(a&b, c)`).
A simple example that showcases the behavior is:
```C
#include <stdio.h>

#ifndef max
    #define max(x, y) ((x > y) ? x : y)
#endif


int main()
{
    int a, b, c;
    a = 1;
    b = 3;
    c = 2;

    // (((1&3) > (2)) ? (1&3) : (2))
    // x = 1&3
    // y = 2
    // (x > y) ? x : y
    // (1&3 > 2) ? 1&3 : 2
    // (1&1) ? 1&3 : 2
    // 1 ?  1&3 : 2

    printf("1&3 = %d\n", a&b);
    printf("max(1&3, 2) = %d\n", max(a&b, c));

    return 0;
}

```

diff_hunk: @@ -3297,7 +3300,9 @@ def check_transformer_n_iter(name, estimator_orig):
         estimator.fit(X, y_)
 
         # These return a n_iter per component.
-        if name in CROSS_DECOMPOSITION:
+        if (name in CROSS_DECOMPOSITION) or (
+            name in ["SVC", "NuSVC"] and len(estimator.classes_) > 2
+        ):
path: sklearn/utils/estimator_checks.py
text: @user1:
```suggestion
        if (name in CROSS_DECOMPOSITION) or (
            # This second condition comes to test the correct
            # exposition of the number of iterations of the
            # libsvm optimization routine.
            name in ["SVC", "NuSVC"] and len(estimator.classes_) > 2
        ):
```

diff_hunk: @@ -3297,7 +3300,9 @@ def check_transformer_n_iter(name, estimator_orig):
         estimator.fit(X, y_)
 
         # These return a n_iter per component.
-        if name in CROSS_DECOMPOSITION:
+        if (name in CROSS_DECOMPOSITION) or (
+            name in ["SVC", "NuSVC"] and len(estimator.classes_) > 2
path: sklearn/utils/estimator_checks.py
text: @user1:
I don't understand why we would get `SVC` and `NuSVC` here. This estimator check (`check_transformer_n_iter`) should only apply to transformers, not classifiers.

It it's not the case, there is problem to fix in the way we routed which checks are applied to which kinds of estimators.

@author:
You are right. We never get `SVC` or `NuSVC` here. I replicated what I did in `check_non_transformer_estimators_n_iter`. Thanks for catching it.

diff_hunk: @@ -3268,9 +3264,16 @@ def check_non_transformer_estimators_n_iter(name, estimator_orig):
 
         set_random_state(estimator, 0)
 
+        if name in ["SVC", "NuSVC"] and estimator.kernel == "precomputed":
+            X = np.dot(X, X.T)
+
         estimator.fit(X, y_)
 
-        assert estimator.n_iter_ >= 1
+        # These return a n_iter per model optimized
+        if name in ["SVC", "NuSVC"] and len(estimator.classes_) > 2:
+            assert np.all(estimator.n_iter_ >= 1)
+        else:
+            assert estimator.n_iter_ >= 1
path: sklearn/utils/estimator_checks.py
text: @user1:
Maybe we could avoid making this test check the name of the estimator special case it and always use:

```python
assert np.all(estimator.n_iter_ >= 1)
```

This should work both for the scalar-valued and array-valued `n_iter_` attributes.