154359Sroberto;#
254359Sroberto;# lr.pl,v 3.1 1993/07/06 01:09:08 jbj Exp
354359Sroberto;#
454359Sroberto;#
554359Sroberto;# Linear Regression Package for perl
654359Sroberto;# to be 'required' from perl
754359Sroberto;#
854359Sroberto;#  Copyright (c) 1992
954359Sroberto;#  Frank Kardel, Rainer Pruy
1054359Sroberto;#  Friedrich-Alexander Universitaet Erlangen-Nuernberg
1154359Sroberto;#
1282498Sroberto;#  Copyright (c) 1997 by
1382498Sroberto;#  Ulrich Windl <Ulrich.Windl@rz.uni-regensburg.de>
1482498Sroberto;#  (Converted to a PERL 5.004 package)
1554359Sroberto;#
1654359Sroberto;#############################################################
1754359Sroberto
1882498Srobertopackage lr;
1982498Sroberto
2054359Sroberto##
2154359Sroberto## y = A + Bx
2254359Sroberto##
2354359Sroberto## B = (n * Sum(xy) - Sum(x) * Sum(y)) / (n * Sum(x^2) - Sum(x)^2)
2454359Sroberto##
2554359Sroberto## A = (Sum(y) - B * Sum(x)) / n
2654359Sroberto##
2754359Sroberto
2854359Sroberto##
2954359Sroberto## interface
3054359Sroberto##
3182498Sroberto;# init(tag);		initialize data set for tag
3282498Sroberto;# sample(x, y, tag);	enter sample
3382498Sroberto;# Y(x, tag);		compute y for given x
3482498Sroberto;# X(y, tag);		compute x for given y
3582498Sroberto;# r(tag);		regression coefficient
3682498Sroberto;# cov(tag);		covariance
3782498Sroberto;# A(tag);
3882498Sroberto;# B(tag);
3982498Sroberto;# sigma(tag);		standard deviation
4082498Sroberto;# mean(tag);
4154359Sroberto#########################
4254359Sroberto
4382498Srobertosub init
4454359Sroberto{
4582498Sroberto    my $self = shift;
4654359Sroberto
4782498Sroberto    $self->{n}   = 0;
4882498Sroberto    $self->{sx}  = 0.0;
4982498Sroberto    $self->{sx2} = 0.0;
5082498Sroberto    $self->{sxy} = 0.0;
5182498Sroberto    $self->{sy}  = 0.0;
5282498Sroberto    $self->{sy2} = 0.0;
5354359Sroberto}
5454359Sroberto
55182007Srobertosub sample($$)
5654359Sroberto{
5782498Sroberto    my $self = shift;
5882498Sroberto    my($_x, $_y) = @_;
5954359Sroberto
6082498Sroberto    ++($self->{n});
6182498Sroberto    $self->{sx}  += $_x;
6282498Sroberto    $self->{sy}  += $_y;
6382498Sroberto    $self->{sxy} += $_x * $_y;
6482498Sroberto    $self->{sx2} += $_x**2;
6582498Sroberto    $self->{sy2} += $_y**2;
6654359Sroberto}
6754359Sroberto
68182007Srobertosub B()
6954359Sroberto{
7082498Sroberto    my $self = shift;
7154359Sroberto
7282498Sroberto    return 1 unless ($self->{n} * $self->{sx2} - $self->{sx}**2);
7382498Sroberto    return ($self->{n} * $self->{sxy} - $self->{sx} * $self->{sy})
7482498Sroberto	/ ($self->{n} * $self->{sx2} - $self->{sx}**2);
7554359Sroberto}
7654359Sroberto
77182007Srobertosub A()
7854359Sroberto{
7982498Sroberto    my $self = shift;
8054359Sroberto
81182007Sroberto    return ($self->{sy} - B() * $self->{sx}) / $self->{n};
8254359Sroberto}
8354359Sroberto
84182007Srobertosub Y()
8554359Sroberto{
8682498Sroberto    my $self = shift;
8754359Sroberto
88182007Sroberto    return A() + B() * $_[$[];
8954359Sroberto}
9054359Sroberto
91182007Srobertosub X()
9254359Sroberto{
9382498Sroberto    my $self = shift;
9454359Sroberto
95182007Sroberto    return ($_[$[] - A()) / B();
9654359Sroberto}
9754359Sroberto
98182007Srobertosub r()
9954359Sroberto{
10082498Sroberto    my $self = shift;
10154359Sroberto
10282498Sroberto    my $s = ($self->{n} * $self->{sx2} - $self->{sx}**2)
10382498Sroberto	  * ($self->{n} * $self->{sy2} - $self->{sy}**2);
10454359Sroberto
10554359Sroberto    return 1 unless $s;
10654359Sroberto
10782498Sroberto    return ($self->{n} * $self->{sxy} - $self->{sx} * $self->{sy}) / sqrt($s);
10854359Sroberto}
10954359Sroberto
110182007Srobertosub cov()
11154359Sroberto{
11282498Sroberto    my $self = shift;
11354359Sroberto
11482498Sroberto    return ($self->{sxy} - $self->{sx} * $self->{sy} / $self->{n})
11582498Sroberto	/ ($self->{n} - 1);
11654359Sroberto}
11754359Sroberto
118182007Srobertosub sigma()
11954359Sroberto{
12082498Sroberto    my $self = shift;
12154359Sroberto
12282498Sroberto    return 0 if $self->{n} <= 1;
12382498Sroberto    return sqrt(($self->{sy2} - ($self->{sy} * $self->{sy}) / $self->{n})
12482498Sroberto		/ ($self->{n}));
12554359Sroberto}
12654359Sroberto
127182007Srobertosub mean()
12854359Sroberto{
12982498Sroberto    my $self = shift;
13054359Sroberto
13182498Sroberto    return 0 if $self->{n} <= 0;
13282498Sroberto    return $self->{sy} / $self->{n};
13354359Sroberto}
13454359Sroberto
13582498Srobertosub new
13682498Sroberto{
13782498Sroberto    my $class = shift;
13882498Sroberto    my $self = {
13982498Sroberto	(n => undef,
14082498Sroberto	 sx => undef,
14182498Sroberto	 sx2 => undef,
14282498Sroberto	 sxy => undef,
14382498Sroberto	 sy => undef,
14482498Sroberto	 sy2 => undef)
14582498Sroberto    };
14682498Sroberto    bless $self, $class;
14782498Sroberto    init($self);
14882498Sroberto    return $self;
14982498Sroberto}
15054359Sroberto
15154359Sroberto1;
152