Regression(3pm)

🦋 NAME 🦋

okay so this is Regression.pm, bestie — the weighted linear regression package, aka the thing that does line + plane fitting ✨

---

💅 SYNOPSIS 💅

first things first, you load it:

use Statistics::Regression;

# Create regression object
my $reg = Statistics::Regression->new( "sample regression", [ "const", "someX", "someY" ] );

# Add data points
$reg->include( 2.0, [ 1.0, 3.0, -1.0 ] );
$reg->include( 1.0, [ 1.0, 5.0, 2.0 ] );
$reg->include( 20.0, [ 1.0, 31.0, 0.0 ] );
$reg->include( 15.0, [ 1.0, 11.0, 2.0 ] );

or, and hear me out, this way:

my %d;
$d{const} = 1.0; $d{someX}= 5.0; $d{someY}= 2.0; $d{ignored}="anything else";
$reg->include( 3.0, \%d );  # names are picked off the Regression specification

🚨 IMPORTANT 🚨 you have to hand it the constant if you want one. it will not do it for you, no cap.

and then, the moment of truth, you print the whole thing:

# Finally, print the result
$reg->print();

which spits out this, literally this:

****************************************************************
Regression 'sample regression'
****************************************************************
Name                         Theta          StdErr     T-stat
[0='const']                 0.2950          6.0512       0.05
[1='someX']                 0.6723          0.3278       2.05
[2='someY']                 1.0688          2.7954       0.38

R^2= 0.808, N= 4
****************************************************************

slay. 🎀

---

okay so the hash input method? lowkey the better one. the way that it works is you can just fill your observation hashes with all your variables, and then run the regression with the exact same code — you only ever change the regression specification in one and only one spot (the new() invocation). you do not need to touch your inputs in the include() statement. like, for real:

my @obs;  ## a global variable.  observations are like: %oneobs= %{$obs[1]};

sub run_regression {
  my $reg = Statistics::Regression->new( $_[0], $_[2] );
  foreach my $obshashptr (@obs) { $reg->include( $_[1], $_[3] ); }
  $reg->print();
}

run_regression("bivariate regression",  $obshashptr->{someY}, [ "const", "someX" ] );
run_regression("trivariate regression",  $obshashptr->{someY}, [ "const", "someX", "someZ" ] );

see?? same sub, two different regressions. iconic behavior tbh 🫶

and obviously, if you wanna do the printing work yourself, the subroutines are right there for you:

my @theta  = $reg->theta();
my @se     = $reg->standarderrors();
my $rsq    = $reg->rsq();
my $adjrsq = $reg->adjrsq();
my $ybar   = $reg->ybar();  ## the average of the y vector
my $sst    = $reg->sst();  ## the sum-squares-total
my $sigmasq= $reg->sigmasq();  ## the variance of the residual
my $k      = $reg->k();   ## the number of variables
my $n      = $reg->n();   ## the number of observations

plus there's some other helper routines hanging around, and a subroutine called linearcombination_variance(). if you don't know what that is — don't use it. that's it, that's the whole warning. 🙈

---

⚠️ BACKGROUND WARNING ⚠️

bestie, you genuinely need to understand OLS regressions to use this package. like, actually. you can get that from an intro college econometrics class and/or from most intermediate college statistics classes. if you do not have this background knowledge, then this package will remain a mystery to you. forever. 🔮

and there is no support for this package — please don't expect any. not a crumb. not a morsel. 💀

---

🌸 DESCRIPTION 🌸

Regression.pm is a multivariate linear regression package. translation: it estimates the c coefficients for a line-fit that looks like this:

y= c(0)*x(0) + c(1)*x1 + c(2)*x2 + ... + c(k)*xk

and you feed it a data set of N observations, each one carrying k independent x variables and one y variable. naturally, N has to be greater than k — and honestly, preferably considerably greater, don't be out here doing the bare minimum. any reasonable undergraduate statistics book will explain what a regression is, so go read one if this is news.

most of the time, you the user will provide a constant ('1') as x(0) for each observation, which lets the regression package fit an intercept. that's the move. 💅

---

🧮 ALGORITHM 🧮

✨ Original Algorithm (ALGOL-60): ✨

W. M. Gentleman, University of Waterloo, "Basic Description For Large, Sparse Or Weighted Linear Least Squares Problems (Algorithm AS 75)," Applied Statistics (1974) Vol 23; No. 3

Gentleman's algorithm is the statistical standard. periodt. insertion of a new observation can be done one observation at any time (WITH A WEIGHT!), and it still only takes a low quadratic time. the storage space requirement is of quadratic order (in the indep variables). a practically infinite number of observations can easily be processed! 🤯

🧠 Internal Data Structures 🧠

R=Rbar is an upperright triangular matrix, kept in normalized form with implicit 1's on the diagonal. D is a diagonal scaling matrix. these correspond to "standard Regression usage" as

                X' X  = R' D R

a back-substitution routine (in thetacov) lets you invert the R matrix — and yes, the inverse is upper-right triangular too, bestie, it's giving symmetry ✨. call this matrix H, that is H=R^(-1).

          (X' X)^(-1) = [(R' D^(1/2)') (D^(1/2) R)]^(-1)
          = [ R^-1 D^(-1/2) ] [ R^-1 D^(-1/2) ]'

---

🐛 BUGS/PROBLEMS 🐛

None known. not a single one. we're clean, we're thriving, we're unbothered. 🧚

😭 Perl Problem 😭

unfortunately, perl is unaware of IEEE number representations. which makes it a whole pain to test whether an observation contains any missing variables (coded as 'NaN' in Regression.pm). it's giving suffering. that's just how it is.

---

📜 VERSION and RECENT CHANGES 📜

• 2007/04/04: Added Coefficient Standard Errors
• 2007/07/01: Added self-test use (if invoked as perl Regression.pm) at the end. cleaned up some print sprintf. changed syntax on new() to eliminate passing K.
• 2007/07/07: allowed passing hash with names to include().

---

👑 AUTHOR 👑

naturally, Gentleman invented this algorithm. it was adaptated by Ivo Welch. Alan Miller (alan\@dmsmelb.mel.dms.CSIRO.AU) pointed out nicer ways to compute the R^2. Ivan Tubert-Brohman helped wrap the module as as a standard CPAN distribution. a whole squad, a whole team, we love to see it 🫶

---

💖 LICENSE 💖

this module is released for free public use under a GPL license.

(C) Ivo Welch, 2001,2004, 2007.

---

and that's the whole page, bestie 🎀 this one's from perl v5.32.1, dated 2021-06-14 — go forth and regress, you've got this 💜