4
$\begingroup$

I am running an experiment several times a day. The experiment has a input parameter that is always greater than 1. The output of the experiment is one of {"OK", "FAIL", "UNKNOWN"}. "UNKNOWN" is very rare and should be ignored.

The output of the experiment is random, but the input parameter influences it. The higher the input parameter is, the less likely the output is to be "OK", but the experiment has a random nature.

My goal is to plot the success rate (# of OK / ( # of OK + # of FAIL)) as a function of the input parameter.

The input parameter is not uniformly distributed, but it is usually close to 2, and very rarely close to 1 or large. I think the x-axis should be transformed by Exp[#] + 1 & to account for this (see partition example below).

Is there a way to plot this from the data in a continuos way? Like the same way that you can interpolate points of data to plot a function from them.

If this is not possible, I was thinking of constructing an array of intervals for the input parameter, then compute (# of OK / ( # of OK + # of FAIL)) for each range and plot that. I have the following code for computing the partition:

{min, max} = {Min @@ #, Max @@ #} & @ data[[All, 1]];
partitionSize = Max[Abs@Log[# - 1] & /@ {min, max}];
PARTS = 20;
partition = 
  Exp[#] + 1 & /@ 
    Range[-partitionSize, partitionSize, 2*partitionSize / parts] // N

If the continuos approach is not possible, would you help me to complete the code so that the data (like the example data below) is grouped in intervals (from the partitition) and then the fraction in each range of the partition is plotted as a bar chart?

Example data:

{{1.49, "OK"}, {1.5, "OK"}, {1.5, "OK"}, {1.52, "OK"}, {1.54, 
  "FAIL"}, {1.54, "FAIL"}, {1.55, "OK"}, {1.57, "OK"}, {1.57, 
  "OK"}, {1.59, "FAIL"}, {1.6, "FAIL"}, {1.6, "UNKNOWN"}, {1.6, 
  "OK"}, {1.61, "FAIL"}, {1.62, "FAIL"}, {1.63, "OK"}, {1.65, 
  "OK"}, {1.65, "FAIL"}, {1.65, "FAIL"}, {1.67, "FAIL"}, {1.68, 
  "OK"}, {1.7, "FAIL"}, {1.7, "OK"}, {1.7, "FAIL"}, {1.7, 
  "OK"}, {1.72, "FAIL"}, {1.72, "OK"}, {1.73, "OK"}, {1.742, 
  "OK"}, {1.75, "FAIL"}, {1.75, "OK"}, {1.75, "FAIL"}, {1.75, 
  "FAIL"}, {1.76, "OK"}, {1.76, "FAIL"}, {1.77, "FAIL"}, {1.78, 
  "OK"}, {1.78, "FAIL"}, {1.79, "OK"}, {1.8, "OK"}, {1.8, 
  "FAIL"}, {1.8, "OK"}, {1.8, "FAIL"}, {1.81, "OK"}, {1.81, 
  "OK"}, {1.83, "FAIL"}, {1.83, "OK"}, {1.83, "OK"}, {1.83, 
  "FAIL"}, {1.83, "OK"}, {1.83, "OK"}, {1.84, "OK"}, {1.84, 
  "OK"}, {1.85, "FAIL"}, {1.85, "FAIL"}, {1.85, "FAIL"}, {1.86, 
  "OK"}, {1.86, "OK"}, {1.87, "FAIL"}, {1.87, "OK"}, {1.87, 
  "OK"}, {1.87, "FAIL"}, {1.88, "OK"}, {1.88, "OK"}, {1.88, 
  "FAIL"}, {1.89, "OK"}, {1.89, "FAIL"}, {1.89, "OK"}, {1.89, 
  "OK"}, {1.89, "OK"}, {1.89, "FAIL"}, {1.89, "OK"}, {1.89, 
  "FAIL"}, {1.89, "FAIL"}, {1.89, "FAIL"}, {1.89, "FAIL"}, {1.89, 
  "FAIL"}, {1.9, "OK"}, {1.9, "OK"}, {1.9, "OK"}, {1.9, 
  "FAIL"}, {1.9, "FAIL"}, {1.9, "OK"}, {1.9, "OK"}, {1.9, 
  "OK"}, {1.9, "OK"}, {1.9, "OK"}, {1.9, "FAIL"}, {1.9, 
  "FAIL"}, {1.9, "FAIL"}, {1.9, "FAIL"}, {1.9, "FAIL"}, {1.9, 
  "FAIL"}, {1.9, "OK"}, {1.9, "OK"}, {1.9, "FAIL"}, {1.9, 
  "OK"}, {1.9, "FAIL"}, {1.9, "OK"}, {1.9, "FAIL"}, {1.9, 
  "FAIL"}, {1.9, "FAIL"}, {1.9, "FAIL"}, {1.9, "OK"}, {1.9, 
  "FAIL"}, {1.9, "FAIL"}, {1.9, "OK"}, {1.9, "OK"}, {1.9, 
  "OK"}, {1.9, "UNKNOWN"}, {1.9, "FAIL"}, {1.9, "FAIL"}, {1.9, 
  "FAIL"}, {1.9, "FAIL"}, {1.9, "OK"}, {1.9, "OK"}, {1.9, 
  "FAIL"}, {1.9, "FAIL"}, {1.9, "FAIL"}, {1.9, "FAIL"}, {1.9, 
  "FAIL"}, {1.9, "FAIL"}, {1.9, "FAIL"}, {1.9, "FAIL"}, {1.9, 
  "OK"}, {1.9, "OK"}, {1.9, "FAIL"}, {1.9, "FAIL"}, {1.9, 
  "FAIL"}, {1.9, "FAIL"}, {1.9, "FAIL"}, {1.9, "FAIL"}, {1.9, 
  "OK"}, {1.9, "FAIL"}, {1.9, "FAIL"}, {1.9, "OK"}, {1.9, 
  "OK"}, {1.9, "OK"}, {1.9095, "OK"}, {1.91, "OK"}, {1.91, 
  "FAIL"}, {1.91, "OK"}, {1.91, "OK"}, {1.91, "OK"}, {1.91, 
  "OK"}, {1.91, "OK"}, {1.91, "OK"}, {1.91, "FAIL"}, {1.91, 
  "OK"}, {1.91, "OK"}, {1.91, "OK"}, {1.91, "FAIL"}, {1.91, 
  "OK"}, {1.91, "FAIL"}, {1.92, "OK"}, {1.94, "OK"}, {1.95, 
  "OK"}, {1.95, "FAIL"}, {1.95, "FAIL"}, {1.95, "OK"}, {1.95, 
  "FAIL"}, {1.95, "OK"}, {1.96, "OK"}, {1.96, "OK"}, {1.96, 
  "FAIL"}, {1.96, "OK"}, {1.96, "FAIL"}, {1.96, "FAIL"}, {1.96, 
  "FAIL"}, {1.96, "FAIL"}, {1.96, "FAIL"}, {1.96, "FAIL"}, {1.96, 
  "OK"}, {1.96, "FAIL"}, {1.96, "OK"}, {1.97, "OK"}, {1.97, 
  "OK"}, {1.97, "OK"}, {1.97, "FAIL"}, {1.97, "FAIL"}, {1.97, 
  "FAIL"}, {1.97, "FAIL"}, {1.97, "FAIL"}, {1.97, "FAIL"}, {1.97065, 
  "OK"}, {1.98, "FAIL"}, {1.98, "FAIL"}, {1.98, "FAIL"}, {1.98, 
  "OK"}, {1.98, "FAIL"}, {1.98, "FAIL"}, {1.98, "OK"}, {1.98, 
  "FAIL"}, {1.98, "OK"}, {1.98, "FAIL"}, {2, "OK"}, {2.05, 
  "OK"}, {2.06, "OK"}, {2.08, "FAIL"}, {2.08, "FAIL"}, {2.14, 
  "FAIL"}, {2.15, "FAIL"}, {2.272, "OK"}, {2.35, "UNKNOWN"}, {2.35, 
  "FAIL"}, {2.35, "OK"}, {2.35, "FAIL"}, {2.38, "FAIL"}, {2.45, 
  "UNKNOWN"}, {2.6, "OK"}, {2.65, "FAIL"}, {2.7, "OK"}, {2.75, 
  "FAIL"}, {2.9, "FAIL"}, {3, "FAIL"}, {3.1, "FAIL"}, {3.21, 
  "FAIL"}, {3.45, "FAIL"}, {3.61, "FAIL"}, {3.61, "FAIL"}, {3.76, 
  "FAIL"}, {3.76, "FAIL"}, {4.32, "FAIL"}, {4.9, "FAIL"}, {5.75, 
  "OK"}, {5.8, "OK"}, {6, "FAIL"}, {6, "FAIL"}, {6.5, 
  "FAIL"}, {6.57, "FAIL"}, {7.3, "FAIL"}, {11.35, "FAIL"}, {12.283, 
  "FAIL"}, {13.14, "FAIL"}, {13.23, "FAIL"}, {21, "FAIL"}}
$\endgroup$
  • $\begingroup$ There is also MinMax[] function which will be handy in your case. {min, max} = MinMax @ data[[All, 1]]; $\endgroup$ – OkkesDulgerci Mar 22 '18 at 16:51
  • $\begingroup$ @OkkesDulgerci thanks for the comment, it certainly is. $\endgroup$ – José D. Mar 22 '18 at 17:38
2
$\begingroup$

This answer will attempt to model the suggested success rate (see question). In what follows, raw is the list of available data.

In order to calculate the success rate (see question) the input data will be grouped into successive bins for the input parameter (see question).

resp = {"OK", "FAIL", "UNKNOWN"};

(* function used to create appropriate dummies  *)
dummy[str_] := Boole@SameQ[#2, str] & @@@ raw /; MemberQ[resp, str];

(* bin the data by x and calc the success rate *)
freq = With[{xs = Part[raw, All, 1], ys = dummy@"OK"},
  {  (* use the mean of x values *)
     Mean[Part[#, All, 1]],
     (* success rate *)
     Plus @@ Part[#, All, -1]/Length[#] 
            (* gather points with abs x difference less than eg 0.05 *)
   } & /@ Gather[Transpose[{xs, ys}], (Abs[#1[[1]] - #2[[1]]] < 0.05 &)]
 ];

freq now holds the binned data along with the success rate.

Before model selection, we are going to consider a couple of additional regressors; plotting the newly acquired data produces the following plot:

Blockquote

It seems that there is considerable concentration of points around a success rate equal to 1 and 0; the additional regressors that will be used will account for those two groups of observations.

data = Join[
   (* success rate = 1, 0 *)
   Boole@Through[{#2 === 1 & @@@ # &, #2 === 0 & @@@ # &}[freq]],
   (* x and success rate (dependent variable) *)
   Transpose[freq]
 ] // Transpose;

data now contains two dummies corresponding to points with success rates equal to 1 and 0 respectively, the input parameter (see question) and the success rate (dependent variable).

There will be several models that will tested; in what follows lmfs contains the result of fitting the different models on their respective data:

lmfs = {
  (* regressing the success rate on the input parameter (a vanilla model) *)
  LinearModelFit[freq, {x}, {x}],

  (* similar to the vanilla model; uses suggestion from the Q *)
  LinearModelFit[freq, {Exp[x]}, {x}],

  (* vanilla model + dummies for extreme cases (success rate = 1\0) *)
  (* note that when d1=d0=0 this is the base case of points with 0 < success rate < 1 *)
  LinearModelFit[data, {d1, d0, x}, {d1, d0, x}],

  (* similar to the previous one with dummies; this allows for cross terms *)
  LinearModelFit[data, {d1, d1 x, d0, d0 x, x}, {d1, d0, x}],

  (* last pair of models: dummies + Exp *)
  LinearModelFit[data, {d1, d0, Exp[x]}, {d1, d0, x}],
  LinearModelFit[data, {d1, d1 Exp[x], d0, d0 Exp[x], Exp[x]}, {d1, d0, x}]
 }

The following table was produced with the code below and it presents the fitted models in ascending order of the AIC criterion:

Join[
  {{"AIC", "RSquared", "AdjustedRSquared", "BasisFunctions"}}, 
  #[{"AIC", "RSquared", "AdjustedRSquared", "BasisFunctions"}] & /@ lmfs // SortBy[#, First] &
 ] // Grid[#, Alignment -> Left] &

Blockquote

It seems that the best fit for the data among the competing models is achieved with the 'hybrid' model with dummies for extreme observations and a transformed (Exp) input parameter (see question).

The following plot presents the (bundled) observations along with the the fitted responses from the 'best' model:

Block[{pts},
 With[{mod = lmfs[[-1]], props = {"Response", "FitResiduals"}, xs = Part[freq, All, 1]},
  pts = Transpose[{xs, #}] & /@ Through[{#[props[[1]]] &, Plus @@ ({1, -1} #[props]) &}[mod]] // N // Chop;

  ListPlot[pts, 
    Frame -> True, 
    FrameLabel -> {x, Row@{success , " ", rate}}, 
    PlotRange -> All, 
    PlotLegends -> {HoldForm@data, fit}]
  ]
 ]

Blockquote

The best fit formula for non-extreme observations (0 < success rate < 1) is given by

(* note how the 'best' model is the last one tried *)
lmfs[[-1]]["BestFit"] /. {d0 -> 0, d1 -> 0}

Blockquote

$\endgroup$
  • $\begingroup$ This is a a great answer, thank you so much. Seems like I still need more data to prove my point. $\endgroup$ – José D. Mar 22 '18 at 21:42
  • $\begingroup$ thank you! do you consider it better than the one you've already accepted? also, could you elaborate more on what the point is? $\endgroup$ – user42582 Mar 22 '18 at 21:45
  • $\begingroup$ I am sorry but I cannot add any interpretation detail. $\endgroup$ – José D. Mar 22 '18 at 21:57
  • $\begingroup$ How are you justifying the regression assumptions of equal variance and (at least roughly) continuous response? $\endgroup$ – JimB Mar 22 '18 at 22:21
  • $\begingroup$ I don't think you need homoskedastic errors to run a linear regression; the dummies are supposed to control for the presence of the extreme values; it's the smooth dependence of the dependent on the regressors that counts, isn't it? $\endgroup$ – user42582 Mar 22 '18 at 23:23
6
$\begingroup$

It sounds like a logistic regression would be one approach to predict the proportion of OK's given the value of your predictor variable x.

(* First get rid of the UNKNOWN's *)
data = Select[data, #[[2]] != "UNKNOWN" &];

(* Turn the OK's to 1's and the FAIL's to 0's *)
data[[All, 2]] = Boole[# == "OK"] & /@ data[[All, 2]];

(* Try a few models *)
lmf0 = LogitModelFit[data, x, x];
lmf1 = LogitModelFit[data, Log[x], x];
lmf2 = LogitModelFit[data, {x, x^2}, x];
lmf3 = LogitModelFit[data, {Log[x], Log[x]^2}, x];

(* Calculate AIC and AIC weight to suggest which model to use *)
(* Here we're simply choosing the model with the smallest AIC value 
   or equivalently the model with the highest AIC weight *)    
aic = {lmf0["AIC"], lmf1["AIC"], lmf2["AIC"], lmf3["AIC"]}
(* {309.32, 307.68, 310.913, 309.317} *)
Δaic = aic - Min[aic]
(* {1.63966, 0., 3.23325, 1.63716} *)
aicWt = Exp[-Δaic/2]/Total[Exp[-Δaic/2]]
(* {0.211769, 0.480738, 0.0954592, 0.212034} *)

(* Plot the results *)
Plot[{lmf0[x], lmf1[x], lmf2[x], lmf3[x]}, {x, Min[data[[All, 1]]], 
  Max[data[[All, 1]]]}, PlotRange -> {All, {0, 0.7}}, 
  PlotStyle -> {{Blue, Thickness[0.001]}, {Red, Thickness[0.01]}, 
  {Green, Thickness[0.001]}, {Orange, Thickness[0.001]}}, 
  PlotLegends -> {"a + b*x", Style["a + b*Log[x]", Bold, 18, Red], 
  "a + b*x + c*x^2", "a + b*Log[x] + c*Log[x]^2"}]

Logistic regression fits

The model named lmf1 has the best fit. The legend shows the logit of the probability of an "OK".

$\endgroup$
  • $\begingroup$ This answer is really helpful, thanks. Is there any other model that does not impose a descending function as answer? I am more interested in what happens for values near 1.9 or so. I am expecting a non-descending function there. $\endgroup$ – José D. Mar 22 '18 at 17:40
  • 1
    $\begingroup$ Logistic regression does not impose a descending function. One can fit other models {Log[x],Log[x]^2,Log[x]^3} for example and/or restrict the predictor values to be within some specified range. $\endgroup$ – JimB Mar 22 '18 at 17:42
4
$\begingroup$

I was thinking of constructing an array of intervals for the input parameter, then compute (# of OK / ( # of OK + # of FAIL)) for each range and plot that.

The measure you describe can be seen as:

  • True Positive Rate (TPR), or

  • Positive Predictive Value (PPV)

from the Receiver Operating Characteristic (ROC) paradigm.

Here is a way to make the plot you describe using the package ROCFunctions.m.

(*Import the package.*)
Import["https://raw.githubusercontent.com/antononcube/MathematicaForPrediction/master/ROCFunctions.m"]

(*Get rid of the UNKNOWN's.*)    
data2 = Select[data, #[[2]] != "UNKNOWN" &];

(*Turn the OKs to 1s and the FAILs to 0s.*)    
data2[[All, 2]] = Boole[# == "OK"] & /@ data2[[All, 2]];

(*Group by parameter.*)
t1 = GroupBy[data2, First -> Last];

(*Simulate a classification process in order to use ROCFunctins.m .*)    
t2 = Map[{{True, True} -> Total[#], {False, False} -> 0, {False, True} -> 0, {True, False} -> Length[#] - Total[#]} &, t1];

(*Add the parameters as ROC parameters.*)    
t3 = KeyValueMap[Association@Join[#2, {"ROCParameter" -> #1}] &, t2];

(*Convert to ROC associations for ROCPlot below.*)    
aROCs = Append[ToROCAssociation[{True, False}, #], "ROCParameter" -> #["ROCParameter"]] & /@ t3;

(*Plot recall vs accuracy.*)
ROCPlot["TPR", "ACC", aROCs, GridLines -> Automatic]

enter image description here

Looking at the comment fragment

[...] I am more interested in what happens for values near 1.9 or so.

in the first answer -- I think the plot below does allow close inspection of what happens around 1.9. (Compare with the plot above.)

rocFuncs = {"TPR", "ACC"};
ListLinePlot[
 Map[Transpose[{Map[#["ROCParameter"] &, aROCs], #}] &, 
  Transpose[Map[Through[ROCFunctions[rocFuncs][#]] &, aROCs]]], 
 Frame -> True, 
 FrameLabel -> Map[Style[#, Larger] &, {"parameter", "rate"}], 
 PlotLegends -> 
  Map[# <> ", " <> (ROCFunctions["FunctionInterpretations"][#]) &, 
   rocFuncs],
 PlotTheme -> "Detailed",
 PlotRange -> {{1.7, 2.1}, {0, 1}}]

enter image description here

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.