LooErrorCSvm.h

Go to the documentation of this file.

/*!
 * 
 *
 * \brief       Leave-one-out error for C-SVMs
 * 
 * 
 *
 * \author      T.Glasmachers
 * \date        2011
 *
 *
 * \par Copyright 1995-2017 Shark Development Team
 * 
 * <BR><HR>
 * This file is part of Shark.
 * <http://shark-ml.org/>
 * 
 * Shark is free software: you can redistribute it and/or modify
 * it under the terms of the GNU Lesser General Public License as published 
 * by the Free Software Foundation, either version 3 of the License, or
 * (at your option) any later version.
 * 
 * Shark is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU Lesser General Public License for more details.
 * 
 * You should have received a copy of the GNU Lesser General Public License
 * along with Shark.  If not, see <http://www.gnu.org/licenses/>.
 *
 */
#ifndef SHARK_OBJECTIVEFUNCTIONS_LOOERRORCSVM_H
#define SHARK_OBJECTIVEFUNCTIONS_LOOERRORCSVM_H


#include <shark/ObjectiveFunctions/AbstractObjectiveFunction.h>
#include <shark/ObjectiveFunctions/Loss/ZeroOneLoss.h>
#include <shark/Algorithms/QP/BoxConstrainedProblems.h>
#include <shark/Algorithms/QP/SvmProblems.h>
#include <shark/LinAlg/CachedMatrix.h>
#include <shark/LinAlg/KernelMatrix.h>
#include <shark/Models/Kernels/KernelExpansion.h>

namespace shark {

///
/// \brief Leave-one-out error, specifically optimized for C-SVMs.
///
template<class InputType, class CacheType = float>
class LooErrorCSvm : public SingleObjectiveFunction
{
public:
    typedef CacheType QpFloatType;
    typedef AbstractKernelFunction<InputType> KernelType;
    typedef LabeledData<InputType, unsigned int> DatasetType;

    /// \brief Constructor.
    LooErrorCSvm(DatasetType const& dataset, KernelType* kernel, bool withOffset)
    : mep_dataset(&dataset)
    , mep_kernel(kernel)
    , m_withOffset(withOffset)
    {
        SHARK_RUNTIME_CHECK(kernel != NULL, "kernel is not allowed to be Null");
        m_features |= HAS_VALUE;
    }


    /// \brief From INameable: return the class name.
    std::string name() const
    { return "LooErrorCSvm"; }
    
    std::size_t numberOfVariables()const{
        return mep_kernel->numberOfParameters()+1;
    }

    /// Evaluate the leave-one-out error for the given parameters. 
    /// These parameters describe the regularization 
    /// constant and the kernel parameters.
    double eval(const RealVector& params){
        QpStoppingCondition stop;
        return eval(params, stop);
    }
    /// Evaluate the leave-one-out error for the given parameters.
    /// These parameters describe the regularization constant and
    /// the kernel parameters.  Furthermore, the stopping
    /// conditions for the solver are passed as an argument.
    double eval(const RealVector& params, QpStoppingCondition &stop){
        this->m_evaluationCounter++;

        double C = params.back();
        mep_kernel->setParameterVector(subrange(params,0,params.size() - 1));
        
        ZeroOneLoss<unsigned int, RealVector> loss;

        // prepare the quadratic program
        typedef KernelMatrix<InputType, QpFloatType> KernelMatrixType;
        typedef CachedMatrix< KernelMatrixType > CachedMatrixType;
        typedef CSVMProblem<CachedMatrixType> SVMProblemType;
        KernelMatrixType km(*mep_kernel, mep_dataset->inputs());
        CachedMatrixType matrix(&km);
        SVMProblemType svmProblem(matrix,mep_dataset->labels(),C);
        std::size_t ell = km.size();
        
        //QpStoppingCondition stop;

        if (m_withOffset)
        {
            // solve the full problem with equality constraint and activated shrinking
            typedef SvmProblem<SVMProblemType> ProblemType;
            ProblemType problem(svmProblem);
            QpSolver< ProblemType > solver(problem);
            solver.solve(stop);
            RealVector alphaFull(problem.dimensions());
            for(std::size_t i = 0; i != problem.dimensions(); ++i){
                alphaFull(i) = problem.alpha(i);
            }
            KernelExpansion<InputType> svm(mep_kernel,mep_dataset->inputs(),true);

            // leave-one-out
            //problem.setShrinking(false);
            double mistakes = 0;
            for (std::size_t i=0; i<ell; i++)
            {
                // use sparseness of the solution:
                if (alphaFull(i) == 0.0) continue;
                problem.deactivateVariable(i);
                
                // solve the reduced problem
                solver.solve(stop);

                // predict using the previously removed example.
                // we need to take into account that the initial problem is solved
                // with shrinking and we thus need to get the initial permutation 
                // for the element index and the unpermuted alpha for the svm
                column(svm.alpha(),0)= problem.getUnpermutedAlpha();
                svm.offset(0) = computeBias(problem);
                std::size_t elementIndex = i;//svmProblem.permutation[i];
                unsigned int target = mep_dataset->element(elementIndex).label;
                mistakes += loss(target, svm(mep_dataset->element(elementIndex).input));
                
                problem.activateVariable(i);
            }
            return mistakes / (double)ell;
        }
        else
        {
            // solve the full problem without equality constraint and activated shrinking
            typedef BoxConstrainedProblem<SVMProblemType> ProblemType;
            ProblemType problem(svmProblem);
            QpSolver< ProblemType > solver(problem);
            solver.solve(stop);
            RealVector alphaFull(problem.dimensions());
            for(std::size_t i = 0; i != problem.dimensions(); ++i){
                alphaFull(i) = problem.alpha(i);
            }
            KernelExpansion<InputType> svm(mep_kernel,mep_dataset->inputs(),false);
            
            // leave-one-out
            //problem.setShrinking(false);
            double mistakes = 0;
            for (std::size_t i=0; i<ell; i++)
            {
                // use sparseness of the solution:
                if (alphaFull(i) == 0.0) continue;
                problem.deactivateVariable(i);
                

                // solve the reduced problem
                solver.solve(stop);

                // predict using the previously removed example.
                // we need to take into account that the initial problem is solved
                // with shrinking and we thus need to get the initial permutation 
                // for the element index and the unpermuted alpha for the svm
                column(svm.alpha(),0)= problem.getUnpermutedAlpha();
                std::size_t elementIndex = i;//svmProblem.permutation[i];
                unsigned int target = mep_dataset->element(elementIndex).label;
                mistakes += loss(target, svm(mep_dataset->element(elementIndex).input));

                problem.activateVariable(i);
            }
            return mistakes / (double)ell;
        }
    }

private:
    /// Compute the SVM offset term (b).
    template<class Problem>
    double computeBias(Problem const& problem){
        double lowerBound = -1e100;
        double upperBound = 1e100;
        double sum = 0.0;
        std::size_t freeVars = 0;
        std::size_t ell = problem.dimensions();
        for (std::size_t i=0; i<ell; i++)
        {
            double value = problem.gradient(i);
            if (problem.alpha(i) == problem.boxMin(i))
            {
                lowerBound = std::max(value,lowerBound);
            }
            else if (problem.alpha(i) == problem.boxMax(i))
            {
                upperBound = std::min(value,upperBound);
            }
            else
            {
                sum += value;
                freeVars++;
            }
        }
        if (freeVars > 0) 
            return sum / freeVars;
        else 
            return 0.5 * (lowerBound + upperBound);
    }
    
    const DatasetType* mep_dataset;
    KernelType* mep_kernel;
    bool m_withOffset;
};


}
#endif