RPN.2.v \ Regressions

\ regressions.rpn -- two variable regressions.
\ Copyright (C) 1999 Eric Laroche.
\ @title "Regressions"
\ @author Eric Laroche
\ @version @(#)$Id: regressions.rpn,v 1.2 1999/05/22 09:14:23 laroche Exp $

\ This program is free software;
\ you can redistribute it and/or modify it.

\ This program is distributed in the hope that it will be useful,
\ but WITHOUT ANY WARRANTY; without even the implied warranty of
\ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.


\ variables:
\ a..f: x, x2, y, y2, xy, n,
\ g..i: lnx, (lnx)2, ylnx,
\ j..l: lny, (lny)2, xlny,
\ m: lnxlny
\ n..o: disable-lnx, disable-lny
\ p: disable-non-linear models
\ q..s: a, b, r2
\ t: parameter validity
\ u: model: best(0),lin,log,exp,pwr
\ v: chosen model


\ clear all variables
\ -->
\ resets all sums to zero
\ resets disable-flags and validity
\ resets models to 'best'
[c] #'21' V{ 0v1+X@ _v} ;

\ check for disable-flags
\ x y --> x y
\ checks x and y sign and p-flag
[f] xn!(0g2r2>!xp|(1Xn)) xo!(0g3r2>!xp|(1Xo)) ;

\ disable parameters validity
\ -->
[d] 0Xt ;

\ calculate parameters
\ x x2 y y2 xy n --> a b r2
\ d(mean((y-yf)**2))/dp = 0
\ a = (x*y-n*xy)/(x**2-n*x2)
\ b = (y-a*x)/n
\ r2 = (n*xy-x*y)**2/((n*x2-x**2)*(n*y2-y**2))
[q] g1g4*g7g7*- g2g1*g8g5*-/ g4g2g4*-g8/ g8r8*g4g7*-g1* g8r6*r5g1*- r7r7*r6g1*-*/ r2r3 ;

\ calculate linear model
\ --> a b r2
[j] xfxexdxcxbxa Cq ;

\ calculate logarithmic model
\ --> a b r2
[g] xn( 0g1g1: xfxixdxcxhxg Cq ) ;

\ calculate exponential model
\ --> a b r2
[e] xo( 0g1g1: xfxlxkxjxbxa Cq ) ;

\ calculate power model
\ --> a b r2
[w] xnxo|( 0g1g1: xfxmxkxjxhxg Cq ) ;

\ call j/g/e/w method
\ i --> a b r2
\ lin(1),log,exp,pwr
[z] c(:Cj:Cg:Ce:Cw) ;

\ display chosen model
\ model -->
[m] c(:D'lin|OK|'d1:D'log|OK|'d1:D'exp|OK|'d1:D'pwr|OK|'d1) ;

\ calculate best model
\ --> a b r2
\ needs to set the model
[k] 0g1 3V{v1+g1Czd1d1 g3g2r2<!(r4r4)d1d1 _v}d1 Xu xuCz xuCm ;

\ unconditionally calculate parameters
\ -->
[u] xvXu xuc(Ck:xuCz) XqXrXs ;

\ calculate parameters
\ -->
[p] xt!( Cu 1Xt ) ;

\ sum single value up or down
\ value down-flag i -->
[t] g1x@ r4r4(-:+) r2X@ ;

\ sum non-logarithmic values up or down
\ down-flag x y --> x y
[n] V Cd g1 v1Ct g1g1* v2Ct g2 v3Ct g2g1* v4Ct g1g3* v5Ct 1 v6Ct ;

\ sum logarithmic values up or down
\ down-flag x y --> x y
[l] V Cf g2xo!(l) g2xn!(l) xn!( g1 v#'7'Ct g1g1* v#'8'Ct g1g5* v#'9'Ct ) xo!( g2 v#'10'Ct g2g1* v#'11'Ct g2g4* v#'12'Ct ) xnxo|!( g1g3* v#'13'Ct ) d1d1 ;

\ sum up or down
\ down-flag x y --> n x y
[s] V vCn vCl xf ;

\ variance
\ sum(v) sum(v**2) n --> var(v)
\ var(v) = (n/(n-1))*mean((v-mean(v))**2)
\ = (sum(v**2)-sum(v)**2/n)/(n-1)
[v] g1*g3/ - r21-/ ;


\ package title
"Regressions"


\ sum up
\ x y --> n x y
"[+]" 0Cs ;

\ sum down
\ x y --> n x y
"[-]" 1Cs ;

\ clear all variables
\ -->
"clr" Cc ;

\ info
\ -->
\ "info" D'Two variable regressions\by Eric Laroche, version 1.2|OK|'d1 ;

~

\ sums
\ --> sum(x) sum(y)
"sum" xc xa ;

\ means
\ --> mean(x) mean(y)
\ mean(v) = sum(v)/n
"mean" xcxf/ xaxf/ ;

\ weighted mean
\ --> wmean(x,y)
\ wmean(v,w) = sum(w*v)/sum(w)
"wm" xexc/ ;

~

\ standard deviations
\ --> sdev(x) sdev(y)
\ sdev(v) = var(v)**(1/2)
"sdev" xfxdxcCvs xfxbxaCvs ;

\ linear modes only flag
\ -->
\ "lin" 1Xp ;

\ all modes flag
\ -->
\ "all" 0Xp ;

\ chose model
\ --> model
"Model" xvUd D'chose model|best|lin|log|exp|pwr|' Xv d1 CdCp xu ;

~

\ x forecast
\ y --> xf
\ xf = (y-b)/a
"x" Cp xuc(:::l:l) xr-xq/ xuc(::e::e) ;

\ y forecast
\ x --> yf
\ yf = a*x+b
"y" Cp xuc(::l::l) xq*xr+ xuc(:::e:e) ;

\ slope
\ --> a
"a" Cp xq ;

\ offset
\ --> b
"b" Cp xr xuc(:::e:e) ;

\ correlation
\ --> r
"r" Cp xss 0xqr2<(n) ;