1/* <module> auc
    2
    3This module computes the Area Under the Receiving Operating Charactersitics and
    4Precision Recall curves using the method of
    5Davis, Jesse, and Mark Goadrich. "The relationship between Precision-Recall
    6and ROC curves."
    7Proceedings of the 23rd international conference on Machine learning. ACM, 2006.
    8
    9@author Fabrizio Riguzzi
   10@license Artistic License 2.0
   11*/
   12
   13:- module(auc,[compute_areas/5,compute_areas_diagrams/5, compute_maxacc/2]).
 compute_areas(+LG:list, -AUCROC:float, -ROC:list, -AUCPR:float, -PR:list) is det
The predicate takes as input 
   32compute_areas(LG,AUCROC,ROC,AUCPR,PR):-
   33  findall(E,member(_- \+(E),LG),Neg),
   34  length(LG,NEx),
   35  length(Neg,NNeg),
   36  NPos is NEx-NNeg,
   37  keysort(LG,LG1),
   38  reverse(LG1,LG2),
   39  compute_pointsroc(LG2,+1e20,0,0,NPos,NNeg,[],ROC),
   40  hull(ROC,0,0,0,AUCROC),
   41  compute_aucpr(LG2,NPos,NNeg,AUCPR,PR).
 compute_areas_diagrams(+LG:list, -AUCROC:float, -ROC:dict, -AUCPR:float, -PR:dict) is det
The predicate takes as input 
   65compute_areas_diagrams(LG,AUCROC,ROC,AUCPR,PR):-
   66  compute_areas(LG,AUCROC,ROC0,AUCPR,PR0),
   67  %  write(ROC0),nl,write(PR0),nl,
   68  ROC = c3{data:_{x:x, rows:[x-'ROC'|ROC0]},
   69    axis:_{x:_{min:0.0,max:1.0,padding:0.0,
   70        tick:_{values:[0.0,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1.0]}},
   71           y:_{min:0.0,max:1.0,padding:_{bottom:0.0,top:0.0},
   72        tick:_{values:[0.0,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1.0]}}}},
   73  PR = c3{data:_{x:x, rows:[x-'PR'|PR0]},
   74    axis:_{x:_{min:0.0,max:1.0,padding:0.0,
   75        tick:_{values:[0.0,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1.0]}},
   76           y:_{min:0.0,max:1.0,padding:_{bottom:0.0,top:0.0},
   77        tick:_{values:[0.0,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1.0]}}}}.
 compute_maxacc(+LG:list, -MaxAcc) is det
The predicate takes as input 
   97compute_maxacc(LG, MaxAcc) :-
   98  findall(E,member(_- \+(E),LG),Neg), %find all the pairs that contain a negative examples
   99  length(LG,NEx),
  100  length(Neg,NNeg),
  101  NPos is NEx-NNeg,
  102  keysort(LG,LG1), % ascending order of the pairs on probabilities
  103  reverse(LG1,LG2), % discending order of the pairs on probabilities
  104  compute_acc_list(LG2, 0, 0, NPos, NNeg, [], AccList),
  105  max_list(AccList, MaxAcc).
 compute_acc_list(+LG:list, +TP:int, +FP:int, +FN:int, +TN:int, +AccList0:list, -AccList:list) is det
The predicate takes as input 
  127compute_acc_list([], TP, FP, FN, TN, AccList0, AccList) :-
  128  Acc is (TP+TN)/(TP+TN+FP+FN),
  129  append(AccList0, [Acc], AccList).
  130
  131compute_acc_list([_P- (\+ _)|T], TP, FP, FN, TN, AccList0, AccList):-!,
  132  Acc is (TP+TN)/(TP+TN+FP+FN),
  133  append(AccList0, [Acc], AccListNew),% append the new accuracy (it creates a new list called AccListNew)
  134  FP1 is FP+1,
  135  TN1 is TN-1,
  136  compute_acc_list(T, TP, FP1, FN, TN1, AccListNew, AccList).
  137
  138compute_acc_list([_P- _|T], TP, FP, FN, TN, AccList0, AccList):-!,
  139  Acc is (TP+TN)/(TP+TN+FP+FN),
  140  append(AccList0, [Acc], AccListNew),
  141  TP1 is TP+1,
  142  FN1 is FN-1,
  143  compute_acc_list(T, TP1, FP, FN1, TN, AccListNew, AccList).
  144
  145
  146
  147
  148compute_pointsroc([],_P0,_TP,_FP,_FN,_TN,P0,P1):-!,
  149  append(P0,[1.0-1.0],P1).
  150
  151compute_pointsroc([P- (\+ _)|T],P0,TP,FP,FN,TN,Po0,Po1):-!,
  152  (P<P0->
  153    FPR is FP/(FP+TN),
  154    TPR is TP/(TP+FN),
  155    append(Po0,[(FPR-TPR)],Po2),
  156    P1=P
  157  ;
  158    Po2=Po0,
  159    P1=P0
  160  ),
  161  FP1 is FP+1,
  162  TN1 is TN-1,
  163  compute_pointsroc(T,P1,TP,FP1,FN,TN1,Po2,Po1).
  164
  165compute_pointsroc([P- _|T],P0,TP,FP,FN,TN,Po0,Po1):-!,
  166  (P<P0->
  167    FPR is FP/(FP+TN),
  168    TPR is TP/(TP+FN),
  169    append(Po0,[FPR-TPR],Po2),
  170    P1=P
  171  ;
  172    Po2=Po0,
  173    P1=P0
  174  ),
  175  TP1 is TP+1,
  176  FN1 is FN-1,
  177  compute_pointsroc(T,P1,TP1,FP,FN1,TN,Po2,Po1).
  178
  179
  180hull([],FPR,TPR,AUC0,AUC1):-
  181  AUC1 is AUC0+(1-FPR)*(1+TPR)/2.
  182
  183
  184hull([FPR1-TPR1|T],FPR,TPR,AUC0,AUC1):-
  185  AUC2 is AUC0+(FPR1-FPR)*(TPR1+TPR)/2,
  186  hull(T,FPR1,TPR1,AUC2,AUC1).
  187
  188compute_aucpr(L,Pos,Neg,A,PR):-
  189  L=[P_0-E|TL],
  190  (E= (\+ _ )->
  191    FP=1,
  192    TP=0,
  193    FN=Pos,
  194    TN is Neg -1
  195  ;
  196    FP=0,
  197    TP=1,
  198    FN is Pos -1,
  199    TN=Neg
  200  ),
  201  compute_curve_points(TL,P_0,TP,FP,FN,TN,Points),
  202  Points=[R0-P0|_TPoints],
  203  (R0=:=0,P0=:=0->
  204    Flag=true
  205  ;
  206    Flag=false
  207  ),
  208  area(Points,Flag,Pos,0,0,0,A,[],PR).
  209
  210compute_curve_points([],_P0,TP,FP,_FN,_TN,[1.0-Prec]):-!,
  211  Prec is TP/(TP+FP).
  212
  213compute_curve_points([P- (\+ _)|T],P0,TP,FP,FN,TN,Pr):-!,
  214  (P<P0->
  215    Prec is TP/(TP+FP),
  216    Rec is TP/(TP+FN),
  217    Pr=[Rec-Prec|Pr1],
  218    P1=P
  219  ;
  220    Pr=Pr1,
  221    P1=P0
  222  ),
  223  FP1 is FP+1,
  224  TN1 is TN-1,
  225  compute_curve_points(T,P1,TP,FP1,FN,TN1,Pr1).
  226
  227compute_curve_points([P- _|T],P0,TP,FP,FN,TN,Pr):-!,
  228  (P<P0->
  229    Prec is TP/(TP+FP),
  230    Rec is TP/(TP+FN),
  231    Pr=[Rec-Prec|Pr1],
  232    P1=P
  233  ;
  234    Pr=Pr1,
  235    P1=P0
  236  ),
  237  TP1 is TP+1,
  238  FN1 is FN-1,
  239  compute_curve_points(T,P1,TP1,FP,FN1,TN,Pr1).
  240
  241area([],_Flag,_Pos,_TPA,_FPA,A,A,PR,PR).
  242
  243area([R0-P0|T],Flag,Pos,TPA,FPA,A0,A,PR0,PR):-
  244 TPB is R0*Pos,
  245  (TPB=:=0->
  246    A1=A0,
  247    FPB=0,
  248    PR2=PR0,
  249    PR=[R0-P0|PR3]
  250  ;
  251    R_1 is TPA/Pos,
  252    (TPA=:=0->
  253      (Flag=false->
  254        P_1=P0,
  255	PR=[0.0-P0|PR3]
  256      ;
  257        P_1=0.0,
  258	PR=[0.0-0.0|PR3]
  259      )
  260    ;
  261      P_1 is TPA/(TPA+FPA),
  262      PR=PR3
  263    ),
  264    FPB is TPB*(1-P0)/P0,
  265    N is TPB-TPA+0.5,
  266    (N<1.0->
  267      append(PR0,[R0-P0],PR2),
  268      A1=A0
  269    ;
  270      interpolate(1,N,Pos,R_1,P_1,TPA,FPA,TPB,FPB,A0,A1,[],PR1),
  271      append(PR0,PR1,PR2)
  272    )
  273  ),
  274  area(T,Flag,Pos,TPB,FPB,A1,A,PR2,PR3).
  275
  276interpolate(I,N,_Pos,_R0,_P0,_TPA,_FPA,_TPB,_FPB,A,A,PR,PR):-I>N,!.
  277
  278interpolate(I,N,Pos,R0,P0,TPA,FPA,TPB,FPB,A0,A,PR0,[R-P|PR]):-
  279  R is (TPA+I)/Pos,
  280  P is (TPA+I)/(TPA+I+FPA+(FPB-FPA)/(TPB-TPA)*I),
  281  A1 is A0+(R-R0)*(P+P0)/2,
  282  I1 is I+1,
  283  interpolate(I1,N,Pos,R,P,TPA,FPA,TPB,FPB,A1,A,PR0,PR)