Skip to content

Commit

Permalink
Add initial cut of Student's T PDF.
Browse files Browse the repository at this point in the history
  • Loading branch information
BrentSeidel committed Jun 13, 2024
1 parent e6cba2b commit b317a81
Show file tree
Hide file tree
Showing 3 changed files with 37 additions and 1 deletion.
17 changes: 17 additions & 0 deletions src/BBS-Numerical-Statistics.adb
Original file line number Diff line number Diff line change
Expand Up @@ -143,4 +143,21 @@ package body BBS.Numerical.Statistics is
end loop;
return 0.0;
end;
--
-- The student's T distrubution is defined as:
-- gamma((nu+1)/2 t*t -(nu+1)/2
-- f(t) = ----------------------- * (1 + ---)
-- sqrt(pi*nu)*gamma(nu/2) nu
--
-- Where nu is the degrees of freedom
--
function studentT_pdf(t : f'Base; nu : Positive) return f'Base is
part1 : f'Base := funct.lngamma2n(nu+1);
part2 : f'Base := elem.log(elem.Sqrt(Ada.Numerics.Pi*f'Base(nu))) +
funct.lngamma2n(nu);
part3 : f'Base := elem.Log(1.0 + (t*t)/f'Base(nu))*(-(f'Base(nu+1)/2.0));
begin
return elem.Exp(part1 + part3 - part2);
end;
--
end;
10 changes: 10 additions & 0 deletions src/BBS-Numerical-Statistics.ads
Original file line number Diff line number Diff line change
Expand Up @@ -44,6 +44,7 @@ package BBS.Numerical.statistics is
--
function chi2_cdf(a, b : F'Base; k, steps : Positive) return F'Base;
--
-- *** Under development - do not use yet. ***
-- Perform a chi^2 test on two sets of data
--
-- This calculates alpha from the expected and observed data sets and then
Expand All @@ -54,6 +55,15 @@ package BBS.Numerical.statistics is
function chi2_test(expected, observed : data_array; k : Positive) return F'Base
with pre => ((expected'First = observed'First) and (expected'Last = observed'Last));
--
-- The student's T distrubution is defined as:
-- gamma((nu+1)/2 t*t -(nu+1)/2
-- f(t) = ----------------------- * (1 + ---)
-- sqrt(pi*nu)*gamma(nu/2) nu
--
-- Where nu is the degrees of freedom
--
function studentT_pdf(t : f'Base; nu : Positive) return f'Base;
--
private
--
-- This is a bit of a hack to get a chi square function that can
Expand Down
11 changes: 10 additions & 1 deletion test/test_stats.adb
Original file line number Diff line number Diff line change
Expand Up @@ -22,6 +22,7 @@ procedure test_stats is
res : linreg.simple_linreg_result;
val : real;
dof : Positive := 50;
sum : real;
begin
Ada.Text_IO.Put_Line("Testing some of the numerical routines.");
res := linreg.simple_linear(data);
Expand Down Expand Up @@ -50,16 +51,24 @@ begin
float_io.Put(stat.chi2_pdf(val, dof), 1, 5, 0);
Ada.Text_IO.Put(" ");
float_io.Put(stat.chi2_cdf(0.0, val, dof, 20), 1, 5, 0);
Ada.Text_IO.Put(" ");
float_io.Put(stat.studentT_pdf(val, dof), 1, 5, 0);
Ada.Text_IO.New_Line;
end loop;
Ada.Text_IO.New_Line;
for n in 0 .. 10 loop
Ada.Text_IO.put(" ");
float_io.put(real(n), 2, 0, 0);
Ada.Text_IO.Put("|");
sum := 0.0;
for k in 0 .. n loop
Ada.Text_IO.Put(" ");
float_io.Put(funct.nchoosek(n, k), 4, 0, 0);
val := funct.nChoosek(n, k);
float_io.Put(val, 4, 0, 0);
sum := sum + val;
end loop;
Ada.Text_IO.Put("| ");
float_io.Put(sum, 4, 0, 0);
Ada.Text_IO.New_Line;
end loop;
end test_stats;

0 comments on commit b317a81

Please sign in to comment.