home
***
CD-ROM
|
disk
|
FTP
|
other
***
search
/
World of Ham Radio 1997
/
WOHR97_AmSoft_(1997-02-01).iso
/
antenna
/
ant_21
/
prog
/
nthreg.bas
< prev
next >
Wrap
BASIC Source File
|
1997-02-01
|
8KB
|
264 lines
REM NTHREG.BAS
cls: PRINT : PRINT " This program finds the coefficients of an Nth order equation using"
PRINT " the method of least squares. From SOME COMMON BASIC PROGRAMS, 3rd"
PRINT " edition, by Lon Poole and Mary Borchers, Osborne/McGraw-Hill,"
PRINT " Berkeley CA, 1979. Authors state original programs were tested on"
PRINT " the Wang 2200 and the Commodore PET (ancient history!). Adapted"
PRINT " where necessary to be compatible with Microsoft QuickBasic 4.5,"
PRINT " modified for more aesthetic and useful output (esp. as an aid in"
PRINT " curve-fitting for programmers), and compiled by Orrin C. Winton WN1Z,"
PRINT " 11/94 and 11/95. The plot routine came from PROGRAMMING WITH BASIC 3rd ed.,"
PRINT " by Byron S. Gottfried, McGraw-Hill, NY, 1986, extensively modified by"
PRINT " O.C. Winton. The menu program is from ADVANCED QuickBASIC by Ken Knecht,"
PRINT " Scott, Foresman and Company, Glenview, Illinois, 1989, and is modified"
PRINT " only slightly.": PRINT
INPUT " <hit any key to continue>", dummy
CLS : PRINT
PRINT " The equation is of the following form:"
PRINT
PRINT " 2 n "
PRINT " y = c + a x + a x + a x "
PRINT " 1 2 n "
PRINT : PRINT
PRINT " where: y = dependent variable "
PRINT " c = constant "
PRINT " a , a ... a = coefficients of independent variables"
PRINT " 1 2 n 2 n "
PRINT " x, x , ... x , respectively "
PRINT
PRINT " The equation coefficients, coefficient of determination, coefficient"
PRINT " of correlation and standard error of estimate are printed."
PRINT
INPUT " <hit any key to continue>", dummy
CLS : PRINT : PRINT " You must provide the x and y coordinates for known data points."
PRINT " Once the equation has been computed you may predict values of y"
PRINT " for given values of x."
PRINT
PRINT " The DIM at line 30 (nthreg.bas) limits the degree of the eqn."
PRINT " You can change this limit according to the following scheme:"
PRINT : PRINT " 30 DIM A(2*D+1), R(D+1, D+2), T(D+2)"
PRINT " where D = maximum degree of equation. E.g.: 30 DIM A(5), R(3,4), T(4)"
PRINT
PRINT
PRINT " Nth-order (or polynomial) regression"
PRINT
REM set limits on degree of equation to A(2D+1), R(D+1,D+2), T(D+2)
REM (where D = maximum degree of equation)
30 DIM a(13), r(7, 8), t(8)
PRINT " degree of equation: ";
INPUT d
PRINT " number of known points: ";
INPUT n
DIM x(30)
DIM y(30)
DIM o(30)
DIM w(30)
a(1) = n
REM enter coordinates of data points
FOR i = 1 TO n
PRINT " x,y of point"; i;
INPUT x(i), y(i)
REM populate matrices with a system of equations
FOR j = 2 TO 2 * d + 1
a(j) = a(j) + x(i) ^ (j - 1)
NEXT j
FOR k = 1 TO d + 1
r(k, d + 2) = t(k) + y(i) * x(i) ^ (k - 1)
t(k) = t(k) + y(i) * x(i) ^ (k - 1)
NEXT k
t(d + 2) = t(d + 2) + y(i) ^ 2
NEXT i
REM solve the system of equations in the matrices
FOR j = 1 TO d + 1
FOR k = 1 TO d + 1
r(j, k) = a(j + k - 1)
NEXT k
NEXT j
FOR j = 1 TO d + 1
k = j
1280 IF r(k, j) <> 0 THEN GOTO 1320
k = k + 1
IF k <= d + 1 THEN GOTO 1280
PRINT " no unique solution"
GOTO 1790
1320 FOR i = 1 TO d + 2
s = r(j, i)
r(j, i) = r(k, i)
r(k, i) = s
NEXT i
z = 1 / r(j, j)
FOR i = 1 TO d + 2
r(j, i) = z * r(j, i)
NEXT i
FOR k = 1 TO d + 1
IF k = j THEN GOTO 1470
z = -r(k, j)
FOR i = 1 TO d + 2
r(k, i) = r(k, i) + z * r(j, i)
NEXT i
1470 NEXT k
NEXT j
PRINT
PRINT " constant = "; r(1, d + 2)
k = r(1, d + 2)
REM print equation coefficients
v = 1
FOR j = 1 TO d
PRINT j; " degree coefficient = "; r(j + 1, d + 2)
u(v) = r(j + 1, d + 2)
v = v + 1
NEXT j
PRINT
REM print as an equation
PRINT
PRINT " f(x) = ";
FOR v = 1 TO d
PRINT "("; u(v); "* x ^"; v; ") + ";
NEXT v
PRINT k
REM compute regression analysis
p = 0
FOR j = 2 TO d + 1
p = p + r(j, d + 2) * (t(j) - a(j) * t(1) / n)
NEXT j
q = t(d + 2) - t(1) ^ 2 / n
z = q - p
i = n - d - 1
PRINT
j = p / q
PRINT " coefficient of determination (r^2) = "; j
PRINT " coefficient of correlation = "; SQR(j)
IF i <= 0 THEN GOTO g
PRINT " standard error of estimate = "; SQR(z / i): GOTO h
g: PRINT " standard error of estimate = "
h: PRINT
REM compute y-coordinate from entered x-coordinate
INPUT " do you wish to do interpolation? enter 'y' for yes"; y$
IF y$ = "y" THEN
PRINT " enter -999 to quit interpolation"
GOTO 1690
ELSE GOTO 1790
END IF
1690 p = r(1, d + 2)
PRINT " x = ";
INPUT x
IF x = -999 THEN GOTO 1790
FOR j = 1 TO d
p = p + r(j + 1, d + 2) * x ^ j
NEXT j
PRINT " y = "; p
PRINT
GOTO 1690
1790 REM
2115 ' ******* Prepare to plot regression equation y-values against
2116 ' user-input x values. Not user-input y-values.
2117 '
2120 FOR i = 1 TO n
2122 o(i) = x(i)
2124 NEXT i
2130 FOR i = 1 TO n
2140 p = r(1, d + 2)
2150 FOR j = 1 TO d
2160 p = p + r(j + 1, d + 2) * o(i) ^ j
2170 NEXT j
2180 w(i) = p
2190 NEXT i
2270 ' ******* Find Largest and Smallest X and Y (the user-input x,y values)
2280 '
2290 XMAX = -100000!: YMAX = -100000!: XMIN = 100000!: YMIN = 100000!
2300 FOR i = 1 TO n
2310 IF x(i) > XMAX THEN XMAX = x(i)
2320 IF y(i) > YMAX THEN YMAX = y(i)
2330 IF x(i) < XMIN THEN XMIN = x(i)
2340 IF y(i) < YMIN THEN YMIN = y(i)
2350 NEXT i
2360
2370 ' ******* Scale the Xs and Ys
2380 '
2390 FOR i = 1 TO n
2400 x(i) = 29 + INT(280 * (x(i) - XMIN) / (XMAX - XMIN))
2410 y(i) = 189 - INT(150 * (y(i) - YMIN) / (YMAX - YMIN))
2420 NEXT i
2470 ' ******* Plot the Graphical Display of user-input x,y values
2480 '
2490 SCREEN 9
2500 LINE (19, 29)-(19, 199): LINE -(319, 199)
2510 FOR IY = 59 TO 179 STEP 20: LINE (19, IY)-(23, IY): NEXT IY
2520 FOR IX = 39 TO 299 STEP 20: LINE (IX, 199)-(IX, 195): NEXT IX
2530 '
2540 FOR i = 1 TO n
2550 PSET (x(i), y(i)): LINE (x(i), y(i) - 2)-(x(i) - 4, y(i) + 2)
2560 LINE -(x(i) + 4, y(i) + 2): LINE -(x(i), y(i) - 2)
2570 NEXT i
2580 '
2590 FOR i = 1 TO n - 1
2600 LINE (x(i), y(i))-(x(i + 1), y(i + 1))
2610 NEXT i
2770 ' ******* Find Largest and Smallest o and w
2780 '
2790 OMAX = -100000!: WMAX = -100000!: OMIN = 100000!: WMIN = 100000!
2800 FOR i = 1 TO n
2810 IF o(i) > OMAX THEN OMAX = o(i)
2820 IF w(i) > WMAX THEN WMAX = w(i)
2830 IF o(i) < OMIN THEN OMIN = o(i)
2840 IF w(i) < WMIN THEN WMIN = w(i)
2850 NEXT i
2860
2870 ' ******* Scale the Os and Ws
2880 '
2890 FOR i = 1 TO n
2900 o(i) = 29 + INT(280 * (o(i) - OMIN) / (OMAX - OMIN))
2910 w(i) = 189 - INT(150 * (w(i) - WMIN) / (WMAX - WMIN))
2920 NEXT i
3000 ' ******* Plot the Regression Equation y-values against
3002 ' user-input x values. Not user-input y-values.
3010 '
3015 COLOR 4
3020 FOR i = 1 TO n - 1
3030 LINE (o(i), w(i))-(o(i + 1), w(i + 1))
3040 NEXT i
3045 COLOR 7
3050 PRINT " f(x) = ";
3060 FOR v = 1 TO d
3070 PRINT "("; u(v); "* x ^"; v; ") + ";
3080 NEXT v
3090 PRINT k
3100 PRINT : PRINT : PRINT : PRINT : PRINT : PRINT : PRINT : PRINT : PRINT : PRINT : PRINT : PRINT : PRINT : PRINT : PRINT : PRINT
3104 PRINT " White curve w/markers: your x,y inputs."
3105 COLOR 4
3106 PRINT " Red curve is eqn's y based on your x-inputs."
3108 PRINT " Eqn's fit is close to perfect if red line overwrites white line."
3109 COLOR 2
3110 INPUT " hit any key to end NTH-ORDER REGRESSION sub-program", dummy
3120 RUN "menu"
3130 END