LazStats: Fix wrong parameter of ProbF() call in BreakDownUnit. Move function ProbF to MathUnit for better testing.

git-svn-id: https://svn.code.sf.net/p/lazarus-ccr/svn@7707 8e941d3f-bd1b-0410-a28a-d453659cc2b4
This commit is contained in:
wp_xxyyzz
2020-09-28 14:05:34 +00:00
parent 4a492958a3
commit 21b70d0447
35 changed files with 221 additions and 232 deletions

View File

@ -1050,7 +1050,7 @@
<Unit127>
<Filename Value="forms\simulations\tprobunit.pas"/>
<IsPartOfProject Value="True"/>
<ComponentName Value="TprobForm"/>
<ComponentName Value="TProbForm"/>
<HasResources Value="True"/>
<ResourceBaseClass Value="Form"/>
<UnitName Value="tProbUnit"/>

View File

@ -91,7 +91,7 @@ var
implementation
uses
Math;
Math, MathUnit;
{ TABCNestedForm }
@ -709,23 +709,23 @@ begin
AReport.Add('SOURCE D.F. SS MS F PROB.');
F := MSA / MSW;
PF := probf(F,dfA,dfwcell);
PF := ProbF(F,dfA,dfwcell);
AReport.Add('A %4D %10.3f%10.3f%10.3f%10.3f', [dfA, SSA, MSA, F, PF]);
F := MSB / MSW;
PF := probf(F,dfBwA,dfwcell);
PF := ProbF(F,dfBwA,dfwcell);
AReport.Add('B(A) %4D %10.3f%10.3f%10.3f%10.3f', [dfBwA, SSB, MSB, F, PF]);
F := MSC / MSW;
PF := probf(F,dfC,dfwcell);
PF := ProbF(F,dfC,dfwcell);
AReport.Add('C %4D %10.3f%10.3f%10.3f%10.3f', [dfC, SSC, MSC, F, PF]);
F := MSAC / MSW;
PF := probf(F,dfAC,dfwcell);
PF := ProbF(F,dfAC,dfwcell);
AReport.Add('AxC %4D %10.3f%10.3f%10.3f%10.3f', [dfAC, SSAC, MSAC, F, PF]);
F := MSBwAC / MSW;
PF := probf(F,dfBwAC,dfwcell);
PF := ProbF(F,dfBwAC,dfwcell);
AReport.Add('B(A)xC %4D %10.3f%10.3f%10.3f%10.3f', [dfBwAC, SSBwAC, MSBwAC, F, PF]);
AReport.Add('w.cells %4D %10.3f%10.3f', [dfwcell, SSW, MSW]);

View File

@ -95,7 +95,7 @@ var
implementation
uses
Math;
Math, MathUnit;
{ TABRAnovaFrm }
@ -559,13 +559,13 @@ begin
FAC := MSAC / MSerrorWithin;
FBC := MSBC / MSerrorWithin;
FABC := MSABC / MSerrorWithin;
ProbA := probf(FA,DFA,DFerrorBetween);
ProbB := probf(FB,DFB,DFerrorBetween);
ProbAB := probf(FAB,DFAB,DFerrorBetween);
ProbC := probf(FC,DFC,DFerrorWithin);
ProbAC := probf(FAC,DFAC,DFerrorWithin);
ProbBC := probf(FBC,DFBC,DFerrorWithin);
ProbABC := probf(FABC,DFABC,DFerrorWithin);
ProbA := ProbF(FA,DFA,DFerrorBetween);
ProbB := ProbF(FB,DFB,DFerrorBetween);
ProbAB := ProbF(FAB,DFAB,DFerrorBetween);
ProbC := ProbF(FC,DFC,DFerrorWithin);
ProbAC := ProbF(FAC,DFAC,DFerrorWithin);
ProbBC := ProbF(FBC,DFBC,DFerrorWithin);
ProbABC := ProbF(FABC,DFABC,DFerrorWithin);
end;
procedure TABRAnovaFrm.Summarize(AReport: TStrings);

View File

@ -120,7 +120,7 @@ var
implementation
uses
Math;
Math, MathUnit;
{ TANCOVAfrm }
@ -780,7 +780,7 @@ begin
MSGroups := SSGroups / df1;
MSError := SSError / df2;
F := MSGroups / MSError;
Prob := probf(F,df1,df2);
Prob := ProbF(F, df1, df2);
AReport.Add('');
AReport.Add('================================================================================');
AReport.Add('');

View File

@ -117,7 +117,7 @@ implementation
uses
Math,
OutputUnit;
MathUnit, OutputUnit;
{ TAxSAnovaFrm }
@ -457,13 +457,13 @@ begin
fd1 := degfree[2];
fd2 := degfree[3];
probf1 := probf(f1, fd1, fd2);
probf1 := ProbF(f1, fd1, fd2);
fd1 := degfree[5];
fd2 := degfree[7];
probf2 := probf(f2, fd1, fd2);
probf2 := ProbF(f2, fd1, fd2);
fd1 := degfree[6];
fd2 := degfree[7];
probf3 := probf(f3, fd1, fd2);
probf3 := ProbF(f3, fd1, fd2);
lReport.Add('Between %5d %10.3f', [degfree[1], ss[1]]);
lReport.Add(' Groups (A) %5d %10.3f %10.3f %10.3f %6.4f', [degfree[2], ss[2], ms[2], f1, probf1]);
lReport.Add(' Subjects w.g.%5d %10.3f %10.3f', [degfree[3], ss[3], ms[3]]);

View File

@ -154,7 +154,8 @@ var
implementation
uses
Math, Utils;
Math,
Utils, MathUnit;
{ TBlksAnovaFrm }
@ -639,7 +640,7 @@ begin
MSDep := SSDep / DFTot;
Omega := (SSF1 - DFF1 * MSErr) / (SSDep + MSErr);
F := MSF1 / MSErr;
ProbF1 := probf(F,DFF1, DFErr);
ProbF1 := ProbF(F,DFF1, DFErr);
MeanDep := MeanDep / N;
end;

View File

@ -78,7 +78,8 @@ var
implementation
uses
Math, Utils;
Math,
Utils, MathUnit;
{ TBNestedAForm }
@ -456,16 +457,16 @@ begin
if RandomBChk.Checked then
begin
F := MSA / MSB;
PF := probf(F, dfA, dfBwA);
PF := ProbF(F, dfA, dfBwA);
end else
begin
F := MSA / MSW;
PF := probf(F, dfA, dfwcell);
PF := ProbF(F, dfA, dfwcell);
end;
AReport.Add('A %4d %10.3f%10.3f%10.3f%10.3f', [dfA, SSA, MSA, F, PF]);
F := MSB / MSW;
PF := probf(F,dfBwA,dfwcell);
PF := ProbF(F,dfBwA,dfwcell);
AReport.Add('B(W) %4d %10.3f%10.3f%10.3f%10.3f', [dfBwA, SSB, MSB, F, PF]);
AReport.Add('w.cells %4d %10.3f%10.3f', [dfwcell, SSW, MSW]);

View File

@ -216,7 +216,7 @@ implementation
uses
Math,
Utils;
Utils, MathUnit;
{ TGLMFrm }
@ -1709,7 +1709,7 @@ begin
F := (R2 / df1) / ((1.0-R2)/ df2)
else
F := 0.0;
FProbF := probf(F,df1,df2);
FProbF := ProbF(F,df1,df2);
AdjR2 := 1.0 - (1.0 - R2) * (totalobs - 1) / df2;
AReport.Add(' R R2 F Prob. > F DF1 DF2 ');

View File

@ -47,7 +47,8 @@ var
implementation
uses
Math, Utils;
Math,
Utils, MathUnit;
const
NO_VALID_NUMBER_ERROR = 'No valid number in row %d of variable "%s"';
@ -298,10 +299,10 @@ begin
fb := MSb / MSwithin;
fc := MSc / MSwithin;
fpartial := MSres / MSwithin;
proba := probf(fa,dfa,dfwithin);
probb := probf(fb,dfb,dfwithin);
probc := probf(fc,dfc,dfwithin);
probpartial := probf(fpartial,dfres,dfwithin);
proba := ProbF(fa,dfa,dfwithin);
probb := ProbF(fb,dfb,dfwithin);
probc := ProbF(fc,dfc,dfwithin);
probpartial := ProbF(fpartial,dfres,dfwithin);
// show ANOVA table results
lReport := TStringList.Create;

View File

@ -131,7 +131,7 @@ var
implementation
uses
Math;
Math, MathUnit;
{ TOneCaseAnovaForm }
@ -655,9 +655,9 @@ begin
FNonAdd := MSNonAdd / MSBalance
else
FNonAdd := 0.0;
ProbF1 := probf(FF1,DFF1,DFErr);
ProbF2 := probf(FF2,DFF2,DFErr);
ProbNonAdd := probf(FNonAdd,1.0,DFBalance);
ProbF1 := ProbF(FF1,DFF1,DFErr);
ProbF2 := ProbF(FF2,DFF2,DFErr);
ProbNonAdd := ProbF(FNonAdd,1.0,DFBalance);
if (ProbF1 > 1.0) then ProbF1 := 1.0;
if (ProbF2 > 1.0) then ProbF2 := 1.0;

View File

@ -86,7 +86,7 @@ var
implementation
uses
Math;
Math, MathUnit;
{ TTtestFrm }
@ -397,11 +397,11 @@ begin
if variance1 > variance2 then
begin
F := variance1 / variance2;
Fp := probf(F,df1,df2);
Fp := ProbF(F, df1, df2);
end else
begin
F := variance2 / variance1;
Fp := probf(F,df2,df1);
Fp := ProbF(F, df2, df1);
end;
lReport.Add('F test for equal variances = %.3f, Probability = %.4f', [F, fp]);
end

View File

@ -84,7 +84,7 @@ var
implementation
uses
Math;
Math, MathUnit;
{ TTwoCorrsFrm }

View File

@ -54,6 +54,9 @@ var
implementation
uses
MathUnit;
{ TWithinANOVAFrm }
procedure TWithinANOVAFrm.ResetBtnClick(Sender: TObject);
@ -237,7 +240,7 @@ begin
ColMeans[i] := ColMeans[i] / count;
end;
f1 := MScols / MSerr; // treatment F statistic
probf1 := probf(f1,dfcols,dferr);
probf1 := ProbF(f1, dfcols, dferr);
// Do reliability terms if requested
if RelChk.Checked then

View File

@ -64,7 +64,7 @@ implementation
uses
Math,
Utils;
Utils, MathUnit;
{ TCannonFrm }
@ -489,18 +489,18 @@ begin
f := (HLTrace * 2.0 * (s * n + 1)) / (s * s * (2.0 * m + s + 1.0));
df1 := s * (2.0 * m + s + 1.0);
df2 := 2.0 * ( s * n + 1.0);
ftestprob := probf(f,df1,df2);
ftestprob := ProbF(f,df1,df2);
lReport.Add('Hotelling-Lawley Trace F-Test %10.4f %2.0f %2.0f %12.4f', [f, df1,df2, ftestprob]);
df2 := s * (2.0 * n + s + 1.0);
f := (Pillia / (s - Pillia)) * ( (2.0 * n + s +1.0) / (2.0 * m + s + 1.0) );
ftestprob := probf(f,df1,df2);
ftestprob := ProbF(f,df1,df2);
lReport.Add('Pillai Trace F-Test %10.4f %2.0f %2.0f %12.4f', [f, df1,df2, ftestprob]);
Roys := Roys * (count - 1 - a_size + b_size)/ a_size ;
df1 := a_size;
df2 := count - 1 - a_size + b_size;
ftestprob := probf(Roys,df1,df2);
ftestprob := ProbF(Roys,df1,df2);
lReport.Add('Roys Largest Root F-Test %10.4f %2.0f %2.0f %12.4f', [Roys, df1, df2, ftestprob]);
lReport.Add('');

View File

@ -68,7 +68,8 @@ var
implementation
uses
Math, Utils;
Math,
Utils, MathUnit;
{ TPartialsFrm }
@ -386,7 +387,7 @@ begin
df1 := TotNoVars - 1 - NoCntrlVars;
df2 := count - TotNoVars;
F := (SemiPart / (1.0 - R2Full)) * (df2 / df1);
Prob := probf(F,df1,df2);
Prob := ProbF(F, df1, df2);
// Report results
lReport.Add('');

View File

@ -63,7 +63,8 @@ var
implementation
uses
Math, Utils;
Math,
Utils, MathUnit;
{ TRMatFrm }
@ -235,7 +236,7 @@ begin
begin
t := Matrix[i-1][j-1] * (sqrt((N-2.0) / (1.0 - (Matrix[i-1][j-1] * Matrix[i-1][j-1]))));
TestMat[i-1,j-1] := t;
Probr := probt(t,N - 2.0);
Probr := ProbT(t,N - 2.0);
TestMat[j-1,i-1] := Probr;
TestMat[i-1,i-1] := 0.0;
end;

View File

@ -59,8 +59,8 @@ type
procedure ANOVA(ListSize: Integer; Freq, Selected, Minimum, Subscript,
Levels, Displace: IntDyneVec; Mean, SS: DblDyneVec; AReport: TStrings);
procedure GetMinMax(AListSize: Integer;
const ASelected: IntDyneVec; var AMinimum, AMaximum: IntDyneVec);
procedure GetLevels(const AMinimum, AMaximum, ALevels, ADisplace: IntDyneVec);
procedure GetMinMax(const AMinimum, AMaximum, ASelected: IntDyneVec);
function Index_Pos(const X, ADisplace: IntDyneVec; AListSize: integer): Integer;
procedure UpdateBtnStates;
@ -79,7 +79,8 @@ implementation
{$R *.lfm}
uses
Math, Utils;
Math,
Utils, MathUnit;
{ TBreakDownFrm }
@ -119,7 +120,6 @@ begin
while true do
begin
//FirstOne:
index := Index_Pos(Subscript, Displace, ListSize);
if Freq[index] > 0 then
begin
@ -127,7 +127,7 @@ begin
for i := 1 to ListSize do
begin
j := Selected[i-1];
AReport.Add('%-10s level = %3d', [
AReport.Add('%-10s level %d', [
OS3MainFrm.DataGrid.Cells[j,0], Minimum[i-1] + subscript[i-1] - 1
]);
end;
@ -155,11 +155,12 @@ begin
MSB := SSB / DF1;
MSW := SSW / DF2;
F := MSB / MSW;
FProb := probf(DF1,DF2,F);
AReport.Add('SOURCE D.F. SS MS F Prob.>F');
AReport.Add('GROUPS %2.0f %8.2f %8.2f %8.3f %6.4f', [DF1, SSB, MSB, F, FProb]);
AReport.Add('WITHIN %2.0f %8.2f %8.2f', [DF2, SSW, MSW]);
AReport.Add('TOTAL %2d %8.2f', [grandsum-1, SST]);
FProb := ProbF(F, DF1, DF2); // wp: was "probf(DF1,DF2,F)";
AReport.Add('SOURCE D.F. SS MS F Prob.>F');
AReport.Add('-------- ---- -------- -------- -------- --------');
AReport.Add('GROUPS %4.0f %8.2f %8.2f %8.3f %8.4f', [DF1, SSB, MSB, F, FProb]);
AReport.Add('WITHIN %4.0f %8.2f %8.2f', [DF2, SSW, MSW]);
AReport.Add('TOTAL %4d %8.2f', [grandsum-1, SST]);
end else
AReport.Add('Insufficient data for ANOVA');
@ -177,7 +178,6 @@ begin
DF1 := 0.0;
DF2 := 0.0;
{ original: }
if ptr1 < 1 then
break;
@ -205,7 +205,7 @@ begin
break;
end;
// do anova for all cells
// Calculate ANOVA for all cells
AReport.Add('ANOVA FOR ALL CELLS');
AReport.Add('');
SST := 0.0;
@ -221,7 +221,7 @@ begin
SST := SST + SS[i];
grandsum := grandsum + Freq[i];
grandsumx := grandsumx + mean[i];
SSW := SSW + (SS[i] - (mean[i] * mean[i] / Freq[i]));
SSW := SSW + (SS[i] - sqr(mean[i]) / Freq[i]);
DF1 := DF1 + 1.0;
DF2 := DF2 + (Freq[i] - 1);
end;
@ -235,11 +235,12 @@ begin
MSB := SSB / DF1;
MSW := SSW / DF2;
F := MSB / MSW;
FProb := probf(DF1, DF2, F);
AReport.Add('SOURCE D.F. SS MS F Prob.>F');
AReport.Add('GROUPS %2.0f %8.2f %8.2f %8.3f %6.4f', [DF1, SSB, MSB, F, FProb]);
AReport.Add('WITHIN %2.0f %8.2f %8.2f', [DF2, SSW, MSW]);
AReport.Add('TOTAL %2d %8.2f', [grandsum-1, SST]);
FProb := ProbF(F, DF1, DF2); // wp: was "probf(DF1, DF2, F)"
AReport.Add('SOURCE D.F. SS MS F Prob.>F');
AReport.Add('-------- ---- -------- -------- -------- --------');
AReport.Add('GROUPS %4.0f %8.2f %8.2f %8.3f %8.4f', [DF1, SSB, MSB, F, FProb]);
AReport.Add('WITHIN %4.0f %8.2f %8.2f', [DF2, SSW, MSW]);
AReport.Add('TOTAL %4d %8.2f', [grandsum-1, SST]);
AReport.Add('');
AReport.Add('FINISHED');
end else
@ -259,9 +260,6 @@ procedure TBreakDownFrm.ComputeBtnClick(Sender: TObject);
label
Label1, Label3, Label4, NextStep;
var
Minimum, Maximum, levels, displace, subscript : IntDyneVec;
Freq : IntDyneVec;
Selected : IntDyneVec;
index, ListSize, length_array : integer;
ptr1, ptr2, sum, grandsum : integer;
xsumtotal, xsqrtotal, grandsumx, grandsumx2, value, SD : double;
@ -270,14 +268,21 @@ var
valstr : string;
dataread : boolean;
selected : IntDyneVec = nil;
freq: IntDyneVec = nil;
minimum: IntDyneVec = nil;
maximum: IntDyneVec = nil;
levels: IntDyneVec = nil;
displace: IntDyneVec = nil;
subscript: IntDyneVec = nil;
mean: DblDyneVec = nil;
variance: DblDyneVec = nil;
stddev: DblDyneVec = nil;
SS: DblDyneVec = nil;
X: Integer;
i, j: integer;
Dependentvar, NoSelected: Integer;
dependentVar, NoSelected: Integer;
tempval: string;
lReport: TStrings;
begin
@ -302,11 +307,6 @@ begin
exit;
end;
// Initialize arrays
SetLength(levels, NoVariables);
SetLength(displace, NoVariables);
SetLength(subscript, NoVariables);
// Get selected variables
SetLength(selected, NoVariables);
for i := 1 to NoSelected do
@ -319,39 +319,19 @@ begin
ListSize := NoSelected;
// Get maximum and minimum levels in each variable
SetLength(minimum, NoVariables);
SetLength(maximum, NoVariables);
GetMinMax(ListSize, Selected, minimum, maximum);
(*
for i := 1 to ListSize do
begin
index := Selected[i-1];
Minimum[i-1] := round(StrToFloat(OS3MainFrm.DataGrid.Cells[index,1]));
Maximum[i-1] := Minimum[i-1];
for j := 1 to NoCases do
begin
if GoodRecord(j,NoSelected,Selected) then
begin
X := round(StrToFloat(OS3MainFrm.DataGrid.Cells[index,j]));
if X < Minimum[i-1] then Minimum[i-1] := X;
if X > Maximum[i-1] then Maximum[i-1] := X;
end;
end;
end;
*)
SetLength(minimum, ListSize);
SetLength(maximum, ListSize);
GetMinMax(minimum, maximum, selected);
// Calculate number of levels for each variable
for i := 1 to ListSize do
levels[i-1] := Maximum[i-1] - Minimum[i-1] + 1;
displace[ListSize-1] := 1;
if ListSize > 1 then
for i := ListSize-1 downto 1 do
displace[i-1] := levels[i] * displace[i];
SetLength(levels, ListSize);
SetLength(displace, ListSize);
GetLevels(minimum, maximum, levels, displace);
// Now, tabulate
length_array := 1;
for i := 1 to ListSize do
length_array := Length_array * levels[i-1];
for i := 0 to ListSize-1 do
length_array := Length_array * levels[i];
// initialize values
SetLength(Freq, length_array+1);
@ -370,10 +350,11 @@ begin
end;
// tabulate
SetLength(subscript, ListSize);
for i := 1 to NoCases do
begin
dataread := false;
if GoodRecord(i,NoSelected,Selected) then
if GoodRecord(i, NoSelected, Selected) then
begin
for j := 1 to ListSize do
begin
@ -386,15 +367,15 @@ begin
end;
if dataread then
begin
j := Index_Pos(subscript,displace,ListSize);
j := Index_Pos(subscript, displace, ListSize);
Freq[j] := Freq[j] + 1;
index := dependentvar;
index := dependentVar;
tempval := Trim(OS3MainFrm.DataGrid.Cells[index,i]);
if tempval <> '' then
begin
value := StrToFloat(tempval);
mean[j] := mean[j] + value;
variance[j] := variance[j] + (value * value);
variance[j] := variance[j] + sqr(value);
end;
end;
end;
@ -605,17 +586,35 @@ begin
end;
procedure TBreakDownFrm.GetMinMax(AListSize: Integer;
const ASelected: IntDyneVec; var AMinimum, AMaximum: IntDyneVec);
procedure TBreakDownFrm.GetLevels(const AMinimum, AMaximum, ALevels, ADisplace: IntDyneVec);
var
i: Integer;
listSize: Integer;
begin
listSize := Length(AMaximum);
for i := 0 to listSize-1 do
ALevels[i] := AMaximum[i] - AMinimum[i] + 1;
ADisplace[listSize-1] := 1;
if listSize > 1 then
for i := listSize-1 downto 1 do
ADisplace[i-1] := ALevels[i] * ADisplace[i];
end;
procedure TBreakDownFrm.GetMinMax(const AMinimum, AMaximum, ASelected: IntDyneVec);
var
i, j, index: Integer;
listSize: Integer;
NoSelected: Integer;
X: Integer;
begin
SetLength(AMinimum, NoVariables);
SetLength(AMaximum, NoVariables);
NoSelected := SelList.Count;
for i := 0 to AListSize-1 do
listSize := Length(AMaximum);
for i := 0 to listSize-1 do
begin
index := ASelected[i];
AMinimum[i] := round(StrToFloat(OS3MainFrm.DataGrid.Cells[index, 1]));
@ -630,8 +629,6 @@ begin
end;
end;
end;
SetLength(AMinimum, AListSize);
SetLengtH(AMaximum, AListSize);
end;
procedure TBreakDownFrm.HelpBtnClick(Sender: TObject);

View File

@ -67,7 +67,8 @@ var
implementation
uses
Math, Utils;
Math,
Utils, MathUnit;
{ TCompareDistFrm }
@ -340,8 +341,8 @@ begin
begin
Xvalue2[i-1] := min2 + (i-1) * incrsize2;
Xvalue2[i] := min2 + (i) * incrsize2;
prob1 := probf(Xvalue2[i-1],df1,df2);
prob2 := probf(Xvalue2[i],df1,df2);
prob1 := ProbF(Xvalue2[i-1],df1,df2);
prob2 := ProbF(Xvalue2[i],df1,df2);
if prob1 > prob2 then
Var2Freq[i-1] := round((prob1-prob2) * Ncases)
else

View File

@ -124,7 +124,8 @@ var
implementation
uses
Math, Utils;
Math,
Utils, MathUnit;
{ TTestScoreFrm }
@ -789,7 +790,7 @@ begin
// R squared values
R2s[i] := 1.0 - (1.0 / CorrMat[i,i]);
W[i] := (R2s[i] / df1) / ((1.0-R2s[i]) / df2);
FProbs[i] := probf(W[i],df1,df2);
FProbs[i] := ProbF(W[i],df1,df2);
valstring := Format('%10s', [ColLabels[i]]);
AReport.Add('%10s %7.3f %7.3f %7.3f %7.3f %4.0f %4.0f', [valstring,sqrt(R2s[i]),R2s[i],W[i],FProbs[i],df1,df2]);

View File

@ -133,7 +133,8 @@ var
implementation
uses
Math, Utils;
Math,
Utils, MathUnit;
{ TBestRegFrm }
@ -593,7 +594,7 @@ begin
ss_res := ( 1.0 - mult_R2) * ss_total ;
ms_res := ss_res / df_res ;
f_test := ms_reg / ms_res ;
prob_f := probf(f_test, df_reg,df_res);
prob_f := ProbF(f_test, df_reg,df_res);
{ Get variance of b coefficients }
AReport.Add('Regression %3d %14.4f %14.4f %14.4f %14.4f', [df_reg, ss_reg, ms_reg, f_test, prob_f]);

View File

@ -82,7 +82,8 @@ var
implementation
uses
Utils, Math;
Math,
Utils, MathUnit;
{ TBlkMregFrm }
@ -318,7 +319,7 @@ begin
df1 := NoIndepVars - pdf1;
df2 := NCases - NoIndepVars - 1;
F := ((R2 - OldR2) / (1.0 - R2)) * df2 / df1;
FProbF := probf(F, df1, df2);
FProbF := ProbF(F, df1, df2);
lReport.Add('F: %26.3f', [F]);
lReport.Add('with probability %10.3f', [FProbF]);
if FProbF < probin then

View File

@ -70,7 +70,7 @@ implementation
uses
Math,
Utils;
Utils, MathUnit;
procedure TLSMregForm.ResetBtnClick(Sender: TObject);
var
@ -298,8 +298,8 @@ begin
df1 := NoIndepVars - pdf1;
df2 := NCases - NoIndepVars - 1;
F := ((R2 - OldR2) / (1.0 - R2)) * df2 / df1;
FProbF := probf(F,df1,df2);
if FProbF < probin then
FProbF := ProbF(F,df1,df2);
if FProbF < probIn then
lReport.Add('Entry requirements met')
else
lReport.Add('Entry requirements not met');

View File

@ -61,7 +61,8 @@ var
implementation
uses
Math, Utils;
Math,
Utils, MathUnit;
{ TSimultFrm }
@ -292,7 +293,7 @@ begin
begin // R squared values
R2s[i-1] := 1.0 - (1.0 / InverseMat[i-1,i-1]);
W[i-1] := (R2s[i-1] / df1) / ((1.0-R2s[i-1]) / df2);
FProbs[i-1] := probf(W[i-1],df1,df2);
FProbs[i-1] := ProbF(W[i-1],df1,df2);
valstring := format('%10s',[ColLabels[i-1]]);
lReport.Add('%10s%10.3f%10.3f%10.3f%10.3f%5.0f%5.0f', [
valstring,sqrt(R2s[i-1]),R2s[i-1],W[i-1],FProbs[i-1],df1,df2

View File

@ -73,7 +73,7 @@ implementation
uses
Math,
Utils;
Utils, MathUnit;
{ TStepFwdFrm }
@ -429,7 +429,7 @@ begin
pdf1 := 1;
pdf2 := NCases - TempNoVars - 1;
PartF := ((NewR2 - R2) * pdf2) / (1.0 - NewR2);
PartProb := probf(PartF, pdf1, pdf2);
PartProb := ProbF(PartF, pdf1, pdf2);
if PartProb < SmallestProb then SmallestProb := PartProb;
if PartProb > LargestProb then LargestProb := PartProb;
lReport.Add('%-10s %6.4f %7.4f %6.4f %3.0f %3.0f', [

View File

@ -107,7 +107,8 @@ var
implementation
uses
Math, Utils;
Math,
Utils, MathUnit;
{ TDiscrimFrm }
@ -500,7 +501,7 @@ begin
GroupMS := GroupSS / DFGroup;
ErrorMS := ErrorSS / DFError;
Fratio := GroupMS / ErrorMS;
prob := probf(Fratio,DFGroup,DFError);
prob := ProbF(Fratio,DFGroup,DFError);
lReport.Add('BETWEEN %4.0f %10.3f %10.3f %10.3f %10.3f', [DFGroup,GroupSS,GroupMS,Fratio,prob]);
lReport.Add('ERROR %4.0f %10.3f %10.3f', [DFError,ErrorSS,ErrorMS]);
lReport.Add('TOTAL %4.0f %10.3f',[DFTotal,TotalSS]);

View File

@ -53,7 +53,7 @@ var
implementation
uses
Math;
Math, MathUnit;
{ TSpearmanFrm }
@ -401,7 +401,7 @@ begin
z := r * sqrt((NCases - 2) / (1.0 - (r * r)));
lReport.Add('t-test value for hypothesis r = 0 is %.3f', [z]);
df := NCases - 2;
Probability := probt(z,df);
Probability := ProbT(z,df);
lReport.Add('Probability > t: %6.4f', [Probability]);
end
else

View File

@ -98,7 +98,7 @@ var
implementation
uses
Math;
Math, MathUnit;
{ TSRHTest }
@ -451,9 +451,9 @@ begin
FF1 := abs(MSF1 / MSErr);
FF2 := abs(MSF2 / MSErr);
FF1F2 := abs(MSF1F2 / MSErr);
ProbF1 := probf(FF1,DFF1,DFErr);
ProbF2 := probf(FF2,DFF2,DFErr);
ProbF1F2 := probf(FF1F2,DFF1F2,DFErr);
ProbF1 := ProbF(FF1,DFF1,DFErr);
ProbF2 := ProbF(FF2,DFF2,DFErr);
ProbF1F2 := ProbF(FF1F2,DFF1F2,DFErr);
if (ProbF1 > 1.0) then ProbF1 := 1.0;
if (ProbF2 > 1.0) then ProbF2 := 1.0;
if (ProbF1F2 > 1.0) then ProbF1F2 := 1.0;

View File

@ -49,7 +49,7 @@ implementation
uses
Math,
Globals, OutputUnit, FunctionsLib;
Globals, OutputUnit, FunctionsLib, MathUnit;
{ TOneSampFrm }

View File

@ -43,7 +43,7 @@ var
implementation
uses
Math;
Math, MathUnit;
{ TFForm }
@ -58,14 +58,12 @@ end;
procedure TFForm.ComputeBtnClick(Sender: TObject);
VAR
F, df1, df2, prob : extended;
outvalue : string;
begin
F := StrToFloat(FEdit.Text);
df1 := StrToFloat(DF1Edit.Text);
df2 := StrToFloat(DF2Edit.Text);
prob := probf(F,df1,df2);
outvalue := format('%6.4f',[prob]);
ProbEdit.Text := outvalue;
prob := ProbF(F,df1,df2);
ProbEdit.Text := Format('%6.4f',[prob]);
end;
procedure TFForm.FormActivate(Sender: TObject);

View File

@ -1,4 +1,4 @@
object TprobForm: TTprobForm
object TProbForm: TTProbForm
Left = 528
Height = 147
Top = 281

View File

@ -5,15 +5,14 @@ unit tProbUnit;
interface
uses
Classes, SysUtils, FileUtil, LResources, Forms, Controls, Graphics, Dialogs,
StdCtrls, ExtCtrls,
Functionslib;
Classes, SysUtils, FileUtil, Forms, Controls, Graphics, Dialogs,
StdCtrls, ExtCtrls;
type
{ TTprobForm }
{ TTProbForm }
TTprobForm = class(TForm)
TTProbForm = class(TForm)
Bevel1: TBevel;
CancelBtn: TButton;
Panel1: TPanel;
@ -36,38 +35,40 @@ type
end;
var
TprobForm: TTprobForm;
TProbForm: TTProbForm;
implementation
{$R *.lfm}
uses
Math;
Math, MathUnit;
{ TTprobForm }
{ TTProbForm }
procedure TTprobForm.ResetBtnClick(Sender: TObject);
procedure TTProbForm.ResetBtnClick(Sender: TObject);
begin
tValueEdit.Text := '';
DFEdit.Text := '';
ProbEdit.Text := '';
end;
procedure TTprobForm.ComputeBtnClick(Sender: TObject);
procedure TTProbForm.ComputeBtnClick(Sender: TObject);
VAR
tvalue, dfvalue, prob : double;
outvalue : string;
begin
tvalue := StrToFloat(tValueEdit.Text);
dfvalue := StrToFloat(DFEdit.Text);
if tvalue >= 0.0 then prob := 0.5 * probt(tvalue,dfvalue);
if tvalue < 0.0 then prob := 1.0 - probt(tvalue,dfvalue) +
if tvalue >= 0.0 then prob := 0.5 * ProbT(tvalue,dfvalue);
if tvalue < 0.0 then prob := 1.0 - ProbT(tvalue,dfvalue) +
(0.5 * probt(tvalue,dfvalue)) ;
if tvalue = 0.0 then prob := 0.50;
outvalue := format('%6.4f',[prob]);
ProbEdit.Text := outvalue;
end;
procedure TTprobForm.FormActivate(Sender: TObject);
procedure TTProbForm.FormActivate(Sender: TObject);
var
w: Integer;
begin
@ -78,8 +79,5 @@ begin
ReturnBtn.Constraints.MinWidth := w;
end;
initialization
{$I tprobunit.lrs}
end.

View File

@ -104,9 +104,7 @@ procedure HomogeneityTest(
implementation
uses
// OutputUnit,
//BlkAnovaUnit,
Utils;
Utils, MathUnit;
procedure Tukey(error_ms : double; { mean squared for residual }
error_df : double; { deg. freedom for residual }
@ -434,7 +432,7 @@ begin
outline := Format('%3d %9.4f ', [i, statistic]);
df1 := 1;
df2 := error_df;
probstat := probf(statistic, round(df1), round(df2)) / 2;
probstat := ProbF(statistic, round(df1), round(df2)) / 2;
outline := outline + Format('%8.3f %5.2f ', [probstat, alpha]);
if probstat < alpha then
outline := outline + 'YES'

View File

@ -10,22 +10,20 @@ uses
MainUnit, dataprocs;
function chisquaredprob(X : double; k : integer) : double;
PROCEDURE matinv(VAR a, vtimesw, v, w: DblDyneMat; n: integer);
FUNCTION sign(a,b: double): double;
FUNCTION isign(a,b : integer): integer;
procedure matinv(VAR a, vtimesw, v, w: DblDyneMat; n: integer);
function sign(a,b: double): double;
function isign(a,b : integer): integer;
function inversez(prob : double) : double;
function zprob(p : double; VAR errorstate : boolean) : double;
function probz(z : double) : double;
function simpsonintegral(a,b : real) : real;
function zdensity(z : real) : real;
function probf(f,df1,df2 : extended) : extended;
FUNCTION alnorm(x : double; upper : boolean): double;
function alnorm(x : double; upper : boolean): double;
procedure ppnd7 (p : double; VAR normal_dev : double; VAR ifault : integer);
function poly(const c: Array of double; nord: integer; x: double): double; // RESULT(fn_val)
procedure swilk (var init : boolean; const x: DblDyneVec; n, n1, n2: integer;
const a: DblDyneVec; var w, pw: double; out ifault: integer);
procedure SVDinverse(VAR a : DblDyneMat; N : integer);
function probt(t,df1 : double) : double;
function inverset(Probt, DF : double) : double;
function inversechi(p : double; k : integer) : double;
function STUDENT(q,v,r : real) : real;
@ -505,51 +503,6 @@ begin
end; (* of normal *)
function probf(f,df1,df2 : extended) : extended;
var
term1, term2, term3, term4, term5, term6 : extended;
FUNCTION gammln(xx: extended): extended;
CONST
stp = 2.50662827465;
// half = 0.5;
one = 1.0;
fpf = 5.5;
VAR
x,tmp,ser: double;
j: integer;
cof: ARRAY [1..6] OF extended;
BEGIN
cof[1] := 76.18009173;
cof[2] := -86.50532033;
cof[3] := 24.01409822;
cof[4] := -1.231739516;
cof[5] := 0.120858003e-2;
cof[6] := -0.536382e-5;
x := xx - 1.0;
tmp := x + fpf;
term1 := ln(tmp);
term2 := (x + 0.5) * term1;
tmp := term2 - tmp;
ser := one;
FOR j := 1 to 6 DO BEGIN
x := x + 1.0;
ser := ser + cof[j] / x
END;
gammln := tmp +ln(stp * ser)
END;
//-----------------------------------------------------------------
begin { fprob function }
if f <= 0.0 then probf := 1.0 else
probf := (betai(0.5*df2,0.5*df1,df2/(df2+df1*f)) +
(1.0-betai(0.5*df1,0.5*df2,df1/(df1+df2/f))))/2.0;
end; // of fprob function
//--------------------------------------------------------------------
{ Algorithm AS66 Applied Statistics (1973) vol.22, no.3
Evaluates the tail area of the standardised normal curve

View File

@ -25,10 +25,12 @@ function BetaI(a,b,x: Double): Extended;
function GammaLn(x: double): Extended;
function tDist(x: Double; N: Integer; OneSided: Boolean): Double;
function tDensity(x: Double; N: Integer): Double;
function tDist(x: Double; DF: Integer; OneSided: Boolean): Double;
function tDensity(x: Double; DF: Integer): Double;
function ProbT(t, DF1: double): double;
function FDensity(x: Double; DF1, DF2: Integer): Double;
function ProbF(F, DF1, DF2: Double): Double;
function Chi2Density(x: Double; N: Integer): Double;
@ -249,23 +251,34 @@ begin
Result := -tmp + ln(2.50662827465 * ser);
end;
// Calculates the (cumulative) t distribution function for N degrees of freedom
function tDist(x: Double; N: Integer; OneSided: Boolean): Double;
// Calculates the (cumulative) t distribution function for DF degrees of freedom
function tDist(x: Double; DF: Integer; OneSided: Boolean): Double;
begin
Result := betai(0.5*N, 0.5, N/(N + sqr(x)));
Result := betai(0.5*DF, 0.5, DF/(DF + sqr(x)));
if OneSided then Result := Result * 0.5;
end;
// Returns the density curve for the t statistic with N degrees of freedom
function tDensity(x: Double; N: Integer): Double;
// Returns the density curve for the t statistic with DF degrees of freedom
function tDensity(x: Double; DF: Integer): Double;
var
factor: Double;
begin
factor := exp(gammaLn((N+1)/2) - gammaLn(N/2)) / sqrt(N * pi);
Result := factor * Power(1 + sqr(x)/N, (1-n)/2);
factor := exp(gammaLn((DF+1)/2) - gammaLn(DF/2)) / sqrt(DF * pi);
Result := factor * Power(1 + sqr(x)/DF, (1-DF)/2);
end;
// Returns the density curve for the F statistic for DF1 and DF2 degrees of freedom.
{ Returns the cumulative probability corresponding to a two-tailed t test with
DF degrees of freedom. }
function ProbT(t, DF1: double): double;
var
F, prob: double;
begin
F := t * t;
Result := ProbF(F, 1.0, DF1);
end;
{ Returns the density curve for the F statistic for DF1 and DF2 degrees of freedom. }
function FDensity(x: Double; DF1, DF2: Integer): Double;
var
ratio1, ratio2, ratio3, ratio4: double;
@ -291,6 +304,23 @@ begin
end;
{ Cumulative probability of the F distribution for given F value and with
DF1 and DF2 degrees of freedom }
function ProbF(F, DF1, DF2: Double): Double;
var
b1, b2: Double;
begin
if f <= 0.0 then
Result := 1.0
else
begin
b1 := betai(0.5*df2, 0.5*df1, df2/(df2+df1*f));
b2 := betai(0.5*df1, 0.5*df2, df1/(df1+df2/f)); // wp here "/f" but in prev line "*f" ????
Result := (b1 + (1.0-b2)) * 0.5;
end;
end;
// Returns the density curve of the chi2 statistic for N degrees of freedom
function Chi2Density(x: Double; N: Integer): Double;
var

View File

@ -142,7 +142,8 @@ procedure matinv(a, vtimesw, v, w: DblDyneMat; n: integer);
implementation
uses
Math, StrUtils, Utils;
Math, StrUtils,
Utils, MathUnit;
procedure GridDotProd(col1, col2: integer; out Product: double; var Ngood: integer);
// Get the cross-product of two vectors
@ -792,7 +793,7 @@ begin
SSreg := SSY - SSres;
R2 := SSreg / SSY;
F := (SSreg / (N - 1)) / (SSres / (NCases - N));
Prob := probf(F,(N-1),(NCases-N));
Prob := ProbF(F, N-1, NCases-N);
AdjR2 := 1.0 - (1.0 - R2) * (NCases - 1) / (NCases - N);
if PrintAll then
begin