testlogitunit.pas

maths lib with source

            XY[I,NF] := 0;
        end
        else
        begin
            
            //
            // label 1
            //
            APVSub(@XY[I][0], 0, NF-1, @C[0], 0, NF-1);
            XY[I,NF] := 1;
        end;
        Inc(I);
    end;
    
    //
    // Train
    //
    MNLTrainH(XY, NP, NF, 2, Info, W, Rep);
    if Info<>1 then
    begin
        ConvErrors := True;
        Exit;
    end;
    
    //
    // Test
    //
    WeakCnt := 0;
    I:=0;
    while I<=NP-1 do
    begin
        APVMove(@X[0], 0, NF-1, @XY[I][0], 0, NF-1);
        MNLProcess(W, X, Y);
        if Y[I mod 2]>0.5 then
        begin
            WeakCnt := WeakCnt+1;
        end;
        J:=0;
        while J<=NC-1 do
        begin
            OtherErrors := OtherErrors or (Y[J]<0) or (Y[J]>1);
            Inc(J);
        end;
        OtherErrors := OtherErrors or (AbsReal(Y[0]+Y[1]-1)>100*MachineEpsilon);
        Inc(I);
    end;
    WeakCorrErrors := WeakCorrErrors or (WeakCnt/NP<2/3);
    
    //
    // Additional tests
    //
    AddTest1(NP, NF, NC, XY, W, OtherErrors);
    OtherErrors := OtherErrors or (MNLClsError(W, XY, NP)<>NP-WeakCnt);
end;


(*************************************************************************
Tests additional functions using already intiialized model.
*************************************************************************)
procedure AddTest1(NP : Integer;
     NF : Integer;
     NC : Integer;
     const XY : TReal2DArray;
     var W : LogitModel;
     var OtherErrors : Boolean);
var
    A : TReal2DArray;
    T : TReal1DArray;
    X : TReal1DArray;
    Y : TReal1DArray;
    W2 : LogitModel;
    Y2 : TReal1DArray;
    I : Integer;
    V : Double;
    S : Double;
    NVars : Integer;
    NClasses : Integer;
begin
    
    //
    // PACK/UNPACK test
    //
    SetLength(T, NC-1+1);
    SetLength(Y, NC-1+1);
    SetLength(Y2, NC-1+1);
    SetLength(X, NF-1+1);
    MNLUnpack(W, A, NVars, NClasses);
    I:=0;
    while I<=NF-1 do
    begin
        X[I] := 2*RandomReal-1;
        Inc(I);
    end;
    S := 0;
    I:=0;
    while I<=NC-2 do
    begin
        V := APVDotProduct(@A[I][0], 0, NF-1, @X[0], 0, NF-1);
        T[I] := Exp(V+A[I,NF]);
        S := S+T[I];
        Inc(I);
    end;
    T[NC-1] := 1;
    S := S+1;
    S := 1/S;
    APVMul(@T[0], 0, NC-1, S);
    MNLProcess(W, X, Y);
    I:=0;
    while I<=NC-1 do
    begin
        OtherErrors := OtherErrors or (AbsReal(Y[I]-T[I])>1000*MachineEpsilon);
        Inc(I);
    end;
    MNLPack(A, NF, NC, W2);
    MNLProcess(W2, X, Y2);
    I:=0;
    while I<=NC-1 do
    begin
        OtherErrors := OtherErrors or (AbsReal(Y2[I]-T[I])>1000*MachineEpsilon);
        Inc(I);
    end;
end;


(*************************************************************************
Tests additional functions.
*************************************************************************)
procedure AddTest2(NF : Integer; NC : Integer; var OtherErrors : Boolean);
var
    A : TReal2DArray;
    XY : TReal2DArray;
    W : LogitModel;
    I : Integer;
    J : Integer;
begin
    SetLength(A, NC-2+1, NF+1);
    I:=0;
    while I<=NC-2 do
    begin
        J:=0;
        while J<=NF do
        begin
            A[I,J] := 0;
            Inc(J);
        end;
        Inc(I);
    end;
    MNLPack(A, NF, NC, W);
    SetLength(XY, NC-1+1, NF+1);
    I:=0;
    while I<=NC-1 do
    begin
        J:=0;
        while J<=NF-1 do
        begin
            XY[I,J] := 2*RandomReal-1;
            Inc(J);
        end;
        XY[I,NF] := I;
        Inc(I);
    end;
    OtherErrors := OtherErrors or (AbsReal(MNLAvgCE(W, XY, NC)-Ln(NC)/Ln(2))>1000*MachineEpsilon);
end;


(*************************************************************************
Tests additional functions.
*************************************************************************)
procedure CVTest(var ConvErrors : Boolean; var OtherErrors : Boolean);
var
    NP : Integer;
    NF : Integer;
    NC : Integer;
    XY : TReal2DArray;
    I : Integer;
    Info : Integer;
    W : LogitModel;
    Rep : MNLReport;
    RelClsExact : Double;
    AvgCEExact : Double;
    NCECVCloser : Integer;
    NThreshold : Integer;
    Pass : Integer;
    PassCount : Integer;
    KSetup : Integer;
begin
    Exit;
end;


(*************************************************************************
Tests error calculation.
*************************************************************************)
procedure ErrTest(var ConvErrors : Boolean; var OtherErrors : Boolean);
var
    NP : Integer;
    NF : Integer;
    NC : Integer;
    XYTrn : TReal2DArray;
    XYTst : TReal2DArray;
    I : Integer;
    J : Integer;
    Info : Integer;
    W : LogitModel;
    Rep : MNLReport;
begin
    
    //
    // Prepare training set such that on test set consisting of 3 examples
    // we have:
    //
    // N    INPUTS      POST.PROB.
    // 0    0 0         (0.25, 0.25, 0.50) for true class 0
    // 1    1 0         (0.25, 0.25, 0.50) for true class 1
    // 2    0 1         (0.25, 0.25, 0.50) for true class 2
    //
    SetLength(XYTrn, 11+1, 2+1);
    I:=0;
    while I<=11 do
    begin
        XYTrn[I,0] := 0;
        XYTrn[I,1] := 0;
        if (I>=4) and (I<8) then
        begin
            XYTrn[I,0] := 1;
        end;
        if (I>=8) and (I<12) then
        begin
            XYTrn[I,1] := 1;
        end;
        XYTrn[I,2] := Min(I mod 4, 2);
        Inc(I);
    end;
    SetLength(XYTst, 2+1, 2+1);
    XYTst[0,0] := 0;
    XYTst[0,1] := 0;
    XYTst[0,2] := 0;
    XYTst[1,0] := 1;
    XYTst[1,1] := 0;
    XYTst[1,2] := 1;
    XYTst[2,0] := 0;
    XYTst[2,1] := 1;
    XYTst[2,2] := 2;
    MNLTrainH(XYTrn, 12, 2, 3, Info, W, Rep);
    if Info<0 then
    begin
        ConvErrors := True;
        Exit;
    end;
    OtherErrors := OtherErrors or (AbsReal(MNLRelClsError(W, XYTst, 3)-2/3)>1000*MachineEpsilon);
    OtherErrors := OtherErrors or (AbsReal(MNLAvgCE(W, XYTst, 3)-(2+2+1)/3)>0.01);
    OtherErrors := OtherErrors or (AbsReal(MNLRMSError(W, XYTst, 3)-Sqrt((2*(Sqr(0.75)+Sqr(0.25)+Sqr(0.5))+Sqr(0.25)+Sqr(0.25)+Sqr(0.50))/9))>0.01);
    OtherErrors := OtherErrors or (AbsReal(MNLAvgError(W, XYTst, 3)-(0.75+0.25+0.50+0.25+0.75+0.50+0.25+0.25+0.50)/9)>0.01);
    OtherErrors := OtherErrors or (AbsReal(MNLAvgRelError(W, XYTst, 3)-(0.75+0.75+0.50)/3)>0.01);
    
    //
    // 3. (1.00, 0.00, 0.00) for true class 0
    // 4. (0.00, 1.00, 0.00) for true class 1
    // 5. (0.00, 0.00, 1.00) for true class 2
    // 6. (0.50, 0.25, 0.25) for true class 0
    // 7. (0.25, 0.50, 0.25) for true class 1
    // 8. (0.25, 0.25, 0.50) for true class 2
    //
end;


(*************************************************************************
Processing functions test
*************************************************************************)
procedure TestProcessing(var ConvErrors : Boolean; var ProcErrors : Boolean);
var
    NVars : Integer;
    NClasses : Integer;
    LM1 : LogitModel;
    LM2 : LogitModel;
    NPoints : Integer;
    XY : TReal2DArray;
    Pass : Integer;
    PassCount : Integer;
    I : Integer;
    J : Integer;
    AllSame : Boolean;
    RLen : Integer;
    Info : Integer;
    Rep : MNLReport;
    X1 : TReal1DArray;
    X2 : TReal1DArray;
    Y1 : TReal1DArray;
    Y2 : TReal1DArray;
    RA : TReal1DArray;
    RA2 : TReal1DArray;
    V : Double;
begin
    PassCount := 100;
    
    //
    // Main cycle
    //
    Pass:=1;
    while Pass<=PassCount do
    begin
        
        //
        // initialize parameters
        //
        NVars := 1+RandomInteger(3);
        NClasses := 2+RandomInteger(3);
        
        //
        // Initialize arrays and data
        //
        NPoints := 10+RandomInteger(50);
        SetLength(X1, NVars-1+1);
        SetLength(X2, NVars-1+1);
        SetLength(Y1, NClasses-1+1);
        SetLength(Y2, NClasses-1+1);
        SetLength(XY, NPoints-1+1, NVars+1);
        I:=0;
        while I<=NPoints-1 do
        begin
            J:=0;
            while J<=NVars-1 do
            begin
                XY[I,J] := 2*RandomReal-1;
                Inc(J);
            end;
            XY[I,NVars] := RandomInteger(NClasses);
            Inc(I);
        end;
        
        //
        // create forest
        //
        MNLTrainH(XY, NPoints, NVars, NClasses, Info, LM1, Rep);
        if Info<=0 then
        begin
            ConvErrors := True;
            Exit;
        end;
        
        //
        // Same inputs leads to same outputs
        //
        I:=0;
        while I<=NVars-1 do
        begin
            X1[I] := 2*RandomReal-1;
            X2[I] := X1[I];
            Inc(I);
        end;
        I:=0;
        while I<=NClasses-1 do
        begin
            Y1[I] := 2*RandomReal-1;
            Y2[I] := 2*RandomReal-1;
            Inc(I);
        end;
        MNLProcess(LM1, X1, Y1);
        MNLProcess(LM1, X2, Y2);
        AllSame := True;
        I:=0;
        while I<=NClasses-1 do
        begin
            AllSame := AllSame and (Y1[I]=Y2[I]);
            Inc(I);
        end;
        ProcErrors := ProcErrors or  not AllSame;
        
        //
        // Same inputs on original forest leads to same outputs
        // on copy created using DFCopy
        //
        UnsetLM(LM2);
        MNLCopy(LM1, LM2);
        I:=0;
        while I<=NVars-1 do
        begin
            X1[I] := 2*RandomReal-1;
            X2[I] := X1[I];
            Inc(I);
        end;
        I:=0;
        while I<=NClasses-1 do
        begin
            Y1[I] := 2*RandomReal-1;
            Y2[I] := 2*RandomReal-1;
            Inc(I);
        end;
        MNLProcess(LM1, X1, Y1);
        MNLProcess(LM2, X2, Y2);
        AllSame := True;
        I:=0;
        while I<=NClasses-1 do
        begin
            AllSame := AllSame and (Y1[I]=Y2[I]);
            Inc(I);
        end;
        ProcErrors := ProcErrors or  not AllSame;
        
        //
        // Same inputs on original forest leads to same outputs
        // on copy created using DFSerialize
        //
        UnsetLM(LM2);
        SetLength(RA, 0+1);
        RA[0] := 0;
        RLen := 0;
        MNLSerialize(LM1, RA, RLen);
        SetLength(RA2, RLen-1+1);
        I:=0;
        while I<=RLen-1 do
        begin
            RA2[I] := RA[I];
            Inc(I);
        end;
        MNLUnserialize(RA2, LM2);
        I:=0;
        while I<=NVars-1 do
        begin
            X1[I] := 2*RandomReal-1;
            X2[I] := X1[I];
            Inc(I);
        end;
        I:=0;
        while I<=NClasses-1 do
        begin
            Y1[I] := 2*RandomReal-1;
            Y2[I] := 2*RandomReal-1;
            Inc(I);
        end;
        MNLProcess(LM1, X1, Y1);
        MNLProcess(LM2, X2, Y2);
        AllSame := True;
        I:=0;
        while I<=NClasses-1 do
        begin
            AllSame := AllSame and (Y1[I]=Y2[I]);
            Inc(I);
        end;
        ProcErrors := ProcErrors or  not AllSame;
        
        //
        // Normalization properties
        //
        I:=0;
        while I<=NVars-1 do
        begin
            X1[I] := 2*RandomReal-1;
            Inc(I);
        end;
        MNLProcess(LM1, X1, Y1);
        V := 0;
        I:=0;
        while I<=NClasses-1 do
        begin
            V := V+Y1[I];
            ProcErrors := ProcErrors or (Y1[I]<0);
            Inc(I);
        end;
        ProcErrors := ProcErrors or (AbsReal(V-1)>1000*MachineEpsilon);
        Inc(Pass);
    end;
end;


(*************************************************************************
Random normal number
*************************************************************************)
function RNormal():Double;
var
    U : Double;
    V : Double;
    S : Double;
    X1 : Double;
    X2 : Double;
begin
    while True do
    begin
        U := 2*RandomReal-1;
        V := 2*RandomReal-1;
        S := Sqr(u)+Sqr(v);
        if (S>0) and (S<1) then
        begin
            S := Sqrt(-2*Ln(S)/S);
            X1 := U*S;
            X2 := V*S;
            Break;
        end;
    end;
    Result := X1;
end;


(*************************************************************************
Random point from sphere
*************************************************************************)
function RSphere(var XY : TReal2DArray; N : Integer; I : Integer):Double;
var
    J : Integer;
    V : Double;
begin
    J:=0;
    while J<=N-1 do
    begin
        XY[I,J] := RNormal;
        Inc(J);
    end;
    V := APVDotProduct(@XY[I][0], 0, N-1, @XY[I][0], 0, N-1);
    V := RandomReal/Sqrt(V);
    APVMul(@XY[I][0], 0, N-1, V);
end;


(*************************************************************************
Unsets model
*************************************************************************)
procedure UnsetLM(var LM : LogitModel);
var
    XY : TReal2DArray;
    Info : Integer;
    Rep : MNLReport;
    I : Integer;
begin
    SetLength(XY, 5+1, 1+1);
    I:=0;
    while I<=5 do
    begin
        XY[I,0] := 0;
        XY[I,1] := 0;
        Inc(I);
    end;
    MNLTrainH(XY, 6, 1, 2, Info, LM, Rep);
    Assert(Info>0);
end;


(*************************************************************************
Silent unit test
*************************************************************************)
function testlogitunit_test_silent():Boolean;
begin
    Result := TestLogit(True);
end;


(*************************************************************************
Unit test
*************************************************************************)
function testlogitunit_test():Boolean;
begin
    Result := TestLogit(False);
end;


end.