*
* First steps in developing the implementation of
* Support Vector Regression in RATS
*
compute nn=41
*
* Gaussian RBF-Kernelfunction
*
function rbfkernel x1 x2 sigma
type real sigma
type real rbfkernel
type vector[real] x1 x2
local MATRIX[REAL] d

comp d = x1-x2
comp rbfkernel = exp(%dot(d,d)*-sigma)
end

* set-up the x and y vectors

declare vector[real] x(nn) y(nn) d(nn)
ewise x(i)=-20.0+(40.0/nn)*(i-1)
ewise y(i)=%if(i==nn/2+1,1.0,sin(x(i))/x(i))+%ran(0.03)
*
* set the parameter for e-insensitive, cost and sigma (Gaussian kernel)
*
comp eps = 0.025
comp cost = 3.0
comp sigma = 0.05

*
* Trying find to estimate the alpha vectors
*
comp nrows = %rows(x)

declare vector[real] alpha(nrows) alphaast(nrows)
declare real costf1 costf2 coatf3 costfunc

nonlin(parmset=pars) alpha alphaast
nonlin(parmset=constraints) %sum(alpha-alphaast)==0.0 alpha>=0.0 alphaast>=0.0 $
      alpha<=cost alphaast<=cost
comp alpha = cost/2*%ones(nrows,1)
comp alphaast = cost/2*%ones(nrows,1)
find(parmset=pars+constraints,method=BFGS,iterations=100,cvcrit=0.1) maximum costfunc
comp costf1 = 0.0
comp costf2 = 0.0
comp costf3 = 0.0
 do i=1, nrows
    do j=1, nrows
       comp x1 = %xrow(X,i)
       comp x2 = %xrow(X,j)
       comp costf1 = costf1 + (alpha(i)-alphaast(i))*(alpha(j)-alphaast(j))*rbfkernel(x1,x2,sigma)
    end
    comp costf2 = costf2 + (alpha(i)+alphaast(i))
    comp costf3 = costf3 + y(i)*(alpha(i)-alphaast(i))
 end
comp costfunc = -0.5*costf1 - eps*costf2 + costf3
disp costfunc -0.5*costf1 -eps*costf2 costf3
end
