Mathematic with Matlab

Linear-Regression

function rFunction = regress( x, y, power )
%REGRESS Fait la régression à la puissance power d'une série de points x,y
    
    % Préparer la matrice
    for i=0:power
        for j=0:power
            exposant=(2*power-i)-j;
            matriceA(i+1,j+1)=sum(x.^exposant)
        end
        matriceR(i+1)=sum(y.*x.^(power-i))
    end
    
    % Résoudre le système
    coeff=transpose(linsolve(matriceA,matriceR(:)))
    exposant=[power:-1:0];
    
    % Imprime les coefficients
    display('Équation:');
    fprintf('\t f(x)=');
    for c=coeff
        if power>0
            fprintf(' %+e*x^%i',c,power);
        else
            fprintf(' %f\n',c);
        end
        power=power-1;
    end

    % Fonction magique,pas important
    rFunction=@(m)arrayfun(@(x) sum(coeff.*(x.^exposant)),m);
    
    % Calcule le coefficient de correllation
    display('R²:');
    fprintf('\t %f \n',(sum((y-mean(y)).^2)-sum((y-rFunction(x)).^2))  /sum((y-mean(y)).^2));
end

Linear-Regression-Golf

function a = golfregress( x, y, p )
    for i=0:p
        for j=0:p
            A(i+1,j+1)=sum(x.^((2*p-i)-j));
        end
        R(i+1)=sum(y.*x.^(p-i));
    end
    a=@(m)arrayfun(@(x)sum((A^-1*R(:))'.*(x.^[p:-1:0])),m);
end

Romberg

function answ = rombergQuad(x,y,p)
    %ROMBERGQUAD Effectue la quadrature de Romberg

    
    % Ceci est une fonction récursive
    % Utiliser avec précotion
    function a=r(k,j)
        if j<=1
            h=(x(length(x))-x(1))/(2^(k-1));
            acc=0;
            m=ceil(length(y)/2^(k-1));
            for i=[2:2^(k-1)]
                acc=acc+y(m*(i-1));
            end
            a=(h/2)*(2*acc+y(1)+y(length(y)));
%            a=(h/2)*(2*y([2:m:y(length(y))])+y(1)+y(length(y)));
        else
           rkim=r(k,j-1);
           a=rkim+(rkim-r(k-1,j-1))/(4^(j-1)-1);
        end
        fprintf('%i,%i : %f\n',k,j,a);
    end
    
    % Effectue les récurtions jusqu'au
    % niveau p
    l=r(p,p);
    answ=l;
end

smartSpline

function answ=cubeSpline(x,y)
    %CUBESPLINE Effectue la spline cubique d'une série de points
    % Utilisation cubeSpline(sérieX, sérieY)

    % Variables
    nPoints=length(x);
    
    % Fonction anonyme pour
    % une écriture plus claire
    h=@(m)  x(m+1)-x(m);
    df=@(m) (y(m+1)-y(m))/h(m);
    dfa=@(m) (y(m)-y(m-1))/h(m-1);
    
    % Construire deux matrices vides
    % A contients les facteurs
    % R contients l'égalitée
    matriceA=zeros(nPoints,nPoints);
    matriceR=zeros(1,nPoints);
    
    % La première et la dernière ligne de la matrice
    % sont simple, puisque C est choisi
    matriceA(1,1)=1;
    matriceA(nPoints,nPoints)=1;
    
    % Construire les autres lignes
    % de la matrice A et R
    for i=2:nPoints-1
        % matrice A
        matriceA(i,i-1)=h(i-1);
        matriceA(i,i)=2*(h(i-1)+h(i));
        matriceA(i,i+1)=h(i);
        
        % matrice R
        matriceR(i)=6*(df(i)-dfa(i));
    end
    % Solutioner la matrice   
    matriceC=linsolve(matriceA,matriceR(:));

    % MATLAB only
    % function rep=lambda(pInter)
    % Both Octave and MATLAB
    function rep=lambda(pInter,x,y,matriceC,nPoints)
        
        % Fonction anonyme pour
        % une écriture plus claire
        % Only Octave need this part
        h=@(m)  x(m+1)-x(m);
        df=@(m) (y(m+1)-y(m))/h(m);
        dfa=@(m) (y(m)-y(m-1))/h(m-1);
        % end of the Octave part
        
        pIndex=1;
        % Evaluer l'interpolation
        % pour chaque points 
        for pEvaluate=pInter
            
            
            % Trouver quel equation
            % il faut employé
            
            % Regarder si le points
            % est dans l'interval
            if pEvaluate>=x(1) && pEvaluate<x(nPoints)
                xIndex=1;
                while ~(pEvaluate>=x(xIndex) && pEvaluate<x(xIndex+1)) && xIndex<nPoints
                    xIndex=xIndex+1;
                end
            elseif pEvaluate<=x(1)
                xIndex=1;
            else
                xIndex=nPoints-1;
            end
           
            % Ajouter le résultats
            C=matriceC;
            
            % Ce n'est pas aussi compliqué que cela le semble
            acc(1)=(C(xIndex)/    (6*h(xIndex)))*(x(xIndex+1)-pInter(pIndex))^3;
            acc(2)=(C(xIndex+1) / (6*h(xIndex)))*(pInter(pIndex)-x(xIndex))^3;
            acc(3)=((y(xIndex)/h(xIndex))-((C(xIndex)*h(xIndex))/6))*(x(xIndex+1)-pInter(pIndex));
            acc(4)=((y(xIndex+1)/h(xIndex))-((C(xIndex+1)*h(xIndex))/6))*(pInter(pIndex)-x(xIndex));
            rep(pIndex)=sum(acc);
            
            pIndex=pIndex+1;
        end
    end

    answ=@(m) lambda(m,x,y,matriceC,nPoints);
end