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 >
BASIC Source File  |  1997-02-01  |  8KB  |  264 lines

  1.      REM   NTHREG.BAS
  2.          cls: PRINT : PRINT " This program finds the coefficients of an Nth order equation using"
  3.      PRINT " the method of least squares.  From SOME COMMON BASIC PROGRAMS, 3rd"
  4.      PRINT " edition, by Lon Poole and Mary Borchers, Osborne/McGraw-Hill,"
  5.      PRINT " Berkeley CA, 1979.  Authors state original programs were tested on"
  6.      PRINT " the Wang 2200 and the Commodore PET (ancient history!).  Adapted"
  7.      PRINT " where necessary to be compatible with Microsoft QuickBasic 4.5,"
  8.      PRINT " modified for more aesthetic and useful output (esp. as an aid in"
  9.      PRINT " curve-fitting for programmers), and compiled by Orrin C. Winton WN1Z,"
  10.      PRINT " 11/94 and 11/95.  The plot routine came from PROGRAMMING WITH BASIC 3rd ed.,"
  11.      PRINT " by Byron S. Gottfried, McGraw-Hill, NY, 1986, extensively modified by"
  12.      PRINT " O.C. Winton.  The menu program is from ADVANCED QuickBASIC by Ken Knecht,"
  13.      PRINT " Scott, Foresman and Company, Glenview, Illinois, 1989, and is modified"
  14.      PRINT " only slightly.": PRINT
  15.      INPUT " <hit any key to continue>", dummy
  16.      CLS : PRINT
  17.      PRINT " The equation is of the following form:"
  18.      PRINT
  19.      PRINT "                             2        n "
  20.      PRINT "      y  =  c  +  a x  +  a x   +  a x  "
  21.      PRINT "                   1       2        n   "
  22.      PRINT : PRINT
  23.      PRINT " where:    y  =  dependent variable       "
  24.      PRINT "           c  =  constant                 "
  25.      PRINT "           a ,  a  ... a   =  coefficients of independent variables"
  26.      PRINT "            1    2      n          2        n                      "
  27.      PRINT "                              x,  x  , ... x  , respectively       "
  28.      PRINT
  29.      PRINT " The equation coefficients, coefficient of determination, coefficient"
  30.      PRINT " of correlation and standard error of estimate are printed."
  31.      PRINT
  32.      INPUT " <hit any key to continue>", dummy
  33.      CLS : PRINT : PRINT " You must provide the x and y coordinates for known data points."
  34.      PRINT " Once the equation has been computed you may predict values of y"
  35.      PRINT " for given values of x."
  36.      PRINT
  37.      PRINT " The DIM at line 30 (nthreg.bas) limits the degree of the eqn."
  38.      PRINT " You can change this limit according to the following scheme:"
  39.      PRINT : PRINT "      30 DIM A(2*D+1), R(D+1, D+2), T(D+2)"
  40.      PRINT " where D  =  maximum degree of equation.  E.g.:  30 DIM A(5), R(3,4), T(4)"
  41.      PRINT
  42.      PRINT
  43.      PRINT " Nth-order (or polynomial) regression"
  44.      PRINT
  45.      REM   set limits on degree of equation to A(2D+1), R(D+1,D+2), T(D+2)
  46.      REM   (where D = maximum degree of equation)
  47. 30       DIM a(13), r(7, 8), t(8)
  48.      PRINT " degree of equation:  ";
  49.      INPUT d
  50.      PRINT " number of known points:  ";
  51.      INPUT n
  52.      DIM x(30)
  53.      DIM y(30)
  54.      DIM o(30)
  55.      DIM w(30)
  56.      a(1) = n
  57.  
  58.      REM   enter coordinates of data points
  59.      FOR i = 1 TO n
  60.          PRINT " x,y of point"; i;
  61.          INPUT x(i), y(i)
  62.  
  63.          REM   populate matrices with a system of equations
  64.          FOR j = 2 TO 2 * d + 1
  65.              a(j) = a(j) + x(i) ^ (j - 1)
  66.          NEXT j
  67.  
  68.          FOR k = 1 TO d + 1
  69.              r(k, d + 2) = t(k) + y(i) * x(i) ^ (k - 1)
  70.              t(k) = t(k) + y(i) * x(i) ^ (k - 1)
  71.          NEXT k
  72.          t(d + 2) = t(d + 2) + y(i) ^ 2
  73.      NEXT i
  74.  
  75.      REM   solve the system of equations in the matrices
  76.      FOR j = 1 TO d + 1
  77.          FOR k = 1 TO d + 1
  78.              r(j, k) = a(j + k - 1)
  79.          NEXT k
  80.      NEXT j
  81.  
  82.      FOR j = 1 TO d + 1
  83.          k = j
  84. 1280     IF r(k, j) <> 0 THEN GOTO 1320
  85.          k = k + 1
  86.          IF k <= d + 1 THEN GOTO 1280
  87.          PRINT " no unique solution"
  88.          GOTO 1790
  89.  
  90. 1320     FOR i = 1 TO d + 2
  91.              s = r(j, i)
  92.              r(j, i) = r(k, i)
  93.              r(k, i) = s
  94.          NEXT i
  95.          z = 1 / r(j, j)
  96.  
  97.          FOR i = 1 TO d + 2
  98.              r(j, i) = z * r(j, i)
  99.          NEXT i
  100.  
  101.          FOR k = 1 TO d + 1
  102.              IF k = j THEN GOTO 1470
  103.              z = -r(k, j)
  104.  
  105.              FOR i = 1 TO d + 2
  106.                  r(k, i) = r(k, i) + z * r(j, i)
  107.              NEXT i
  108. 1470     NEXT k
  109.      NEXT j
  110.  
  111.      PRINT
  112.      PRINT "             constant  =  "; r(1, d + 2)
  113.      k = r(1, d + 2)
  114.  
  115.      REM   print equation coefficients
  116.      v = 1
  117.      FOR j = 1 TO d
  118.          PRINT j; " degree coefficient  =  "; r(j + 1, d + 2)
  119.          u(v) = r(j + 1, d + 2)
  120.          v = v + 1
  121.      NEXT j
  122.      PRINT
  123.  
  124.      REM   print as an equation
  125.      PRINT
  126.      PRINT " f(x) = ";
  127.      FOR v = 1 TO d
  128.          PRINT "("; u(v); "* x ^"; v; ") + ";
  129.      NEXT v
  130.      PRINT k
  131.  
  132.  
  133.      REM   compute regression analysis
  134.      p = 0
  135.      FOR j = 2 TO d + 1
  136.          p = p + r(j, d + 2) * (t(j) - a(j) * t(1) / n)
  137.      NEXT j
  138.  
  139.      q = t(d + 2) - t(1) ^ 2 / n
  140.      z = q - p
  141.      i = n - d - 1
  142.      PRINT
  143.      j = p / q
  144.      PRINT " coefficient of determination (r^2) = "; j
  145.      PRINT " coefficient of correlation         = "; SQR(j)
  146.      IF i <= 0 THEN GOTO g
  147.      PRINT " standard error of estimate         = "; SQR(z / i): GOTO h
  148. g:       PRINT " standard error of estimate         = "
  149. h:       PRINT
  150.  
  151.      REM   compute y-coordinate from entered x-coordinate
  152.      INPUT " do you wish to do interpolation?  enter 'y' for yes"; y$
  153.      IF y$ = "y" THEN
  154.          PRINT " enter -999 to quit interpolation"
  155.          GOTO 1690
  156.      ELSE GOTO 1790
  157.      END IF
  158.  
  159. 1690     p = r(1, d + 2)
  160.      PRINT " x = ";
  161.      INPUT x
  162.      IF x = -999 THEN GOTO 1790
  163.      FOR j = 1 TO d
  164.          p = p + r(j + 1, d + 2) * x ^ j
  165.      NEXT j
  166.      PRINT " y = "; p
  167.      PRINT
  168.      GOTO 1690
  169. 1790     REM
  170.  
  171.  
  172. 2115 ' ******* Prepare to plot regression equation y-values against
  173. 2116 '         user-input x values.  Not user-input y-values.
  174. 2117 '
  175. 2120 FOR i = 1 TO n
  176. 2122     o(i) = x(i)
  177. 2124 NEXT i
  178.  
  179.  
  180. 2130 FOR i = 1 TO n
  181. 2140     p = r(1, d + 2)
  182. 2150     FOR j = 1 TO d
  183. 2160         p = p + r(j + 1, d + 2) * o(i) ^ j
  184. 2170     NEXT j
  185. 2180     w(i) = p
  186. 2190 NEXT i
  187.  
  188.  
  189. 2270 ' ******* Find Largest and Smallest X and Y (the user-input x,y values)
  190. 2280 '
  191. 2290 XMAX = -100000!: YMAX = -100000!: XMIN = 100000!: YMIN = 100000!
  192. 2300 FOR i = 1 TO n
  193. 2310    IF x(i) > XMAX THEN XMAX = x(i)
  194. 2320    IF y(i) > YMAX THEN YMAX = y(i)
  195. 2330    IF x(i) < XMIN THEN XMIN = x(i)
  196. 2340    IF y(i) < YMIN THEN YMIN = y(i)
  197. 2350 NEXT i
  198. 2360
  199. 2370 ' ******* Scale the Xs and Ys
  200. 2380 '
  201. 2390 FOR i = 1 TO n
  202. 2400     x(i) = 29 + INT(280 * (x(i) - XMIN) / (XMAX - XMIN))
  203. 2410     y(i) = 189 - INT(150 * (y(i) - YMIN) / (YMAX - YMIN))
  204. 2420 NEXT i
  205.  
  206. 2470 ' ******* Plot the Graphical Display of user-input x,y values
  207. 2480 '
  208. 2490 SCREEN 9
  209. 2500 LINE (19, 29)-(19, 199): LINE -(319, 199)
  210. 2510 FOR IY = 59 TO 179 STEP 20: LINE (19, IY)-(23, IY): NEXT IY
  211. 2520 FOR IX = 39 TO 299 STEP 20: LINE (IX, 199)-(IX, 195): NEXT IX
  212. 2530 '
  213. 2540 FOR i = 1 TO n
  214. 2550     PSET (x(i), y(i)): LINE (x(i), y(i) - 2)-(x(i) - 4, y(i) + 2)
  215. 2560     LINE -(x(i) + 4, y(i) + 2): LINE -(x(i), y(i) - 2)
  216. 2570 NEXT i
  217. 2580 '
  218. 2590 FOR i = 1 TO n - 1
  219. 2600     LINE (x(i), y(i))-(x(i + 1), y(i + 1))
  220. 2610 NEXT i
  221.  
  222. 2770 ' ******* Find Largest and Smallest o and w
  223. 2780 '
  224. 2790 OMAX = -100000!: WMAX = -100000!: OMIN = 100000!: WMIN = 100000!
  225. 2800 FOR i = 1 TO n
  226. 2810    IF o(i) > OMAX THEN OMAX = o(i)
  227. 2820    IF w(i) > WMAX THEN WMAX = w(i)
  228. 2830    IF o(i) < OMIN THEN OMIN = o(i)
  229. 2840    IF w(i) < WMIN THEN WMIN = w(i)
  230. 2850 NEXT i
  231. 2860
  232. 2870 ' ******* Scale the Os and Ws
  233. 2880 '
  234. 2890 FOR i = 1 TO n
  235. 2900     o(i) = 29 + INT(280 * (o(i) - OMIN) / (OMAX - OMIN))
  236. 2910     w(i) = 189 - INT(150 * (w(i) - WMIN) / (WMAX - WMIN))
  237. 2920 NEXT i
  238.  
  239. 3000 ' ******* Plot the Regression Equation y-values against
  240. 3002 '         user-input x values.  Not user-input y-values.
  241. 3010 '
  242. 3015 COLOR 4
  243. 3020 FOR i = 1 TO n - 1
  244. 3030     LINE (o(i), w(i))-(o(i + 1), w(i + 1))
  245. 3040 NEXT i
  246.  
  247. 3045 COLOR 7
  248. 3050 PRINT " f(x) = ";
  249. 3060 FOR v = 1 TO d
  250. 3070    PRINT "("; u(v); "* x ^"; v; ") + ";
  251. 3080 NEXT v
  252. 3090 PRINT k
  253.  
  254. 3100 PRINT : PRINT : PRINT : PRINT : PRINT : PRINT : PRINT : PRINT : PRINT : PRINT : PRINT : PRINT : PRINT : PRINT : PRINT : PRINT
  255. 3104 PRINT " White curve w/markers: your x,y inputs."
  256. 3105 COLOR 4
  257. 3106 PRINT " Red curve is eqn's y based on your x-inputs."
  258. 3108 PRINT " Eqn's fit is close to perfect if red line overwrites white line."
  259. 3109 COLOR 2
  260. 3110 INPUT " hit any key to end NTH-ORDER REGRESSION sub-program", dummy
  261. 3120 RUN "menu"
  262. 3130 END
  263.  
  264.