home *** CD-ROM | disk | FTP | other *** search
/ ftp.barnyard.co.uk / 2015.02.ftp.barnyard.co.uk.tar / ftp.barnyard.co.uk / cpm / walnut-creek-CDROM / CPM / TURBOPAS / MAPSTATF.LBR / DESCRPT.PZS / DESCRPT.ÐAS
Text File  |  2000-06-30  |  9KB  |  259 lines
  1. (* Multivariate Analysis Package - Descriptive Statistics Module
  2.    Copyright 1985 Douglas L. Anderton.  This program may be freely
  3.    circulated so long as it is not sold for profit and any charge does
  4.    not exceed costs of reproduction *)
  5.  
  6. Program Descrpt(Input,Output);
  7. Label initer, ngroup;
  8. Const
  9.    N=100;
  10. Type
  11.    SUBS = 1..N;
  12.    RVEC = Array [SUBS] Of Real;
  13.    IVEC = Array [SUBS] Of Integer;
  14.    S8 = Array [SUBS] Of String[8];
  15.    S = String[8];
  16. Var
  17.    dfile, ofile: Text;
  18.    sel : IVEC;
  19.    t1, t2, t3, vr, gold, wht, mu, stdv, ster, vrnc, skw, krt, rng: Real;
  20.    nc, miss, vars, mm1, mm2, mm3, mm4, vmin, vmax : RVEC;
  21.    varn : S8;
  22.    grp, fg, i, j, k, nv, dv, ot, gp, wt, wgt : Integer;
  23.    ngp : Boolean;
  24.    yn : S;
  25.  
  26. Procedure openfiles(Var dfile, ofile:Text; Var ot:Integer);
  27. Var
  28.    dfl, ofl:String[12];
  29. Begin
  30.    ClrScr; Writeln(' *** DESCRPT: DESCRIPTIVE STATISTICS ***');
  31.    Writeln; Write('Name of the data file? ');
  32.    Readln(dfl); Assign(dfile,dfl); Reset(dfile);
  33.    Write('Name of the output file? ');
  34.    Readln(ofl); Assign(ofile,ofl); Rewrite(ofile);
  35.    ot:=0;
  36.    If(ofl='CON:') Or (ofl='con:') Then ot:=1;
  37.    If(ofl='LST:') Or (ofl='lst:') Then ot:=2;
  38.    If(ot=2) Then
  39.       Begin
  40.       Writeln(ofile,'Multivariate Analysis Package (1.6) - ',
  41.          'Program: DESCRPT, Datafile: ',dfl); Writeln(ofile);
  42.       End;
  43.    End; (* Of openfiles *)
  44.  
  45. Procedure wait(ot:Integer);
  46. Begin
  47.   If ot=1 Then Begin
  48.     Write('- Press any key to continue -'); While Not KeyPressed Do; ClrScr;
  49.     End;
  50.   End; (* of wait *)
  51.  
  52. Procedure selvar(Var sel:IVEC; Var varn:S8; Var miss:RVEC; 
  53.                  Var gp, wt, dv, nv:Integer);
  54. Var
  55.    cfile:Text;
  56.    cfl:String[12];
  57.    i,j,f:Integer;
  58.    mis:Real;
  59.    van:String[8];
  60. Begin
  61.    Write('Name of the codebook file (or NONE)? '); Readln(cfl);
  62.    If (cfl<>'NONE') And (cfl<>'none') Then f:=1 Else f:=0;
  63.    If f=1 Then Begin Assign(cfile,cfl); Reset(cfile); End;
  64.    Writeln;
  65.    Write('How many variables in data file? '); Readln(nv);
  66.    Write('Number of variables to use in DESCRPT? '); Readln(dv);
  67.    For i:=1 To dv Do
  68.      Begin
  69.      Write('Column number for variable ',i,'? '); Readln(sel[i]);
  70.      Str(sel[i]:3,varn[i]); miss[i]:=-1E37;
  71.      End;
  72.    Write('Of these Column numbers which is weight (0=none)? '); Readln(wt);
  73.    Write('Of these Column numbers which is grouping (0=none)? '); Readln(gp);
  74.    If gp<>0 Then Begin
  75.      Writeln('Note: data assumed sorted on grouping variable');
  76.      Writeln('      histograms not allowed with grouped data');
  77.      End;
  78.    If f=1 Then Begin
  79.      For j:=1 to nv Do Begin
  80.        mis:=-1E37;
  81.        Readln(cfile,f,van,mis);
  82.        For i:=1 to dv Do
  83.          If f=sel[i] Then Begin
  84.            varn[i]:=van; miss[i]:=mis;
  85.            Writeln('Col: ',sel[i],'  Name: ',varn[i],' Missing: ',miss[i]:6);
  86.            End;     
  87.        End;
  88.      Close(cfile);
  89.      wait(1);
  90.      End;
  91.    End; (* Of selvar *)
  92.  
  93. Procedure getcase(Var vars:RVEC; sel:IVEC; nv, dv:Integer; Var dfile:Text);
  94. Var
  95.    i, j:Integer;
  96.    x:Real;
  97. Begin
  98.    For i:=1 To nv Do
  99.       Begin
  100.       Read(dfile,x);
  101.       For j:=1 To dv Do If(sel[j]=i) Then vars[j]:=x;
  102.       End;
  103.    End; (* Of getcase *)
  104.  
  105. Procedure histogram(mx, mn, iv:Real; vnum:Integer; varn:S; miss:Real;
  106.                     Var dfile, ofile:Text);
  107. Var
  108.    l, k, tot:Integer;
  109.    x, p, c:Real;
  110.    s:String[1];
  111.    freq:Array [1..200] Of Integer;
  112. Begin
  113.    (* reread data For frequencies *)
  114.    Reset(dfile);
  115.    tot:=0;
  116.    FillChar(freq,2*200,0);
  117.    While Not Eof(dfile) Do
  118.       Begin
  119.       For l:=1 To vnum Do Read(dfile,x);
  120.       Readln(dfile);
  121.       If(x <> miss) Then
  122.         Begin
  123.         l:=Trunc(((x-mn)/((mx-mn)/iv))+1);
  124.         If l>Trunc(iv) Then l:=Trunc(iv);
  125.         freq[l]:=freq[l]+1; tot:=tot+1;
  126.         End;
  127.       End;
  128.    (* now print graphics *)
  129.    s:='#'; c:=0.0;
  130.    ClrScr; Writeln(ofile);
  131.    Writeln(ofile, 'Histogram for ', varn,' with ',Round(iv),' intervals');
  132.    Writeln(ofile);
  133.    Writeln(ofile,'From     To       Freq  Pct    Cum   --10--20--30--40--50',
  134.       '--60--70--80--90-100');
  135.    mx:=(mx-mn)/iv;
  136.    For l:=1 To Trunc(iv) Do
  137.      Begin
  138.      p:=(freq[l]/tot)*100.0;
  139.      c:=c+p;
  140.      Write(ofile,mn+((l-1)*mx):9:3,mn+(l*mx):9:3,
  141.         freq[l]:4,p:7:2,c:7:2,'|');
  142.      For k:=1 To Round(0.4*p) Do Write(ofile,s);
  143.      Writeln(ofile);
  144.      End;
  145.    Writeln(ofile,'                                     --10--20--30--40--50',
  146.       '--60--70--80--90-100');
  147.    End; (* Of histogram *)
  148.  
  149. Begin (* main *)
  150.    openfiles(dfile, ofile, ot);
  151.    selvar(sel, varn, miss, gp, wt, dv, nv);
  152.    (* intialize *)
  153.    fg:=0;
  154.    wgt:=0; grp:=0;
  155.    For i:=1 To dv Do
  156.       Begin
  157.       If wt=sel[i] Then wgt:=i;
  158.       If gp=sel[i] Then grp:=i;
  159.       End;
  160. initer:
  161.    FillChar(nc,6*N,0); FillChar(mm1,6*N,0); FillChar(mm2,6*N,0);
  162.    FillChar(mm3,6*N,0); FillChar(mm4,6*N,0); 
  163.    (* accumulate *)
  164.    Writeln;
  165.    k:=0;
  166.    While Not EOF(dfile) Do
  167.       Begin
  168.       k:=k+1;
  169.       If fg=0 Then getcase(vars, sel, nv, dv, dfile);
  170.       fg:=0;
  171.       If Frac(k/10)=0.0 Then Write('+');
  172.       ngp:=False;
  173.       If Not EOF(dfile) Then
  174.          Begin
  175.          If((grp<>0) And (k>1)) Then ngp:=Not(vars[grp]=gold);
  176.          gold:=vars[grp];
  177.          If ngp Then Goto ngroup;
  178.          If(((wgt<>0) And (vars[wgt]<>miss[wgt])) Or (wgt=0)) Then
  179.          For j:=1 To dv Do
  180.             Begin
  181.             If vars[j]<>miss[j] Then
  182.                Begin
  183.                wht:=1;
  184.                If(wgt<>0) And (j<>wgt) Then wht:=vars[wgt];
  185.                If nc[j]=0.0 Then
  186.                   Begin
  187.                   vmax[j]:=vars[j]; vmin[j]:=vars[j];
  188.                   End;
  189.                If vars[j]>vmax[j] Then vmax[j]:=vars[j];
  190.                If vars[j]<vmin[j] Then vmin[j]:=vars[j];
  191.                (* Spicer one-pass algorithm *)
  192.                t1:=nc[j];
  193.                nc[j]:=nc[j]+wht;
  194.                vr:=(vars[j]-mm1[j])*wht/nc[j];
  195.                t2:=Sqr(vr);
  196.                t3:=Sqr(wht);
  197.                mm4[j]:=mm4[j]-vr*mm3[j]*4.0+t2*mm2[j]*6.0 +
  198.                          ((Sqr(nc[j])-3*wht*t1)/(t3*wht))*
  199.                          Sqr(t2)*nc[j]*t1;
  200.                mm3[j]:=mm3[j]-vr*mm2[j]*3.0+(nc[j]*t1/
  201.                          t3)*(nc[j]-(wht+wht))*t2*vr;
  202.                mm2[j]:=mm2[j]+(nc[j]*t1/wht)*t2;
  203.                mm1[j]:=mm1[j]+vr;
  204.                End;
  205.             End;
  206.          End;
  207.       End;
  208. ngroup:
  209.    (* compute And Output stats*)
  210.    For j:=1 To dv Do
  211.       Begin
  212.       ClrScr;
  213.       Writeln(ofile,'Descriptive Statistics for ',varn[j]); Writeln(ofile);
  214.       If nc[j] > 1 Then
  215.         Begin
  216.         mu:=mm1[j]; vrnc:=mm2[j]/(nc[j]-1.0);
  217.         stdv:=Sqrt(vrnc); ster:=stdv/Sqrt(nc[j]);
  218.         rng:=vmax[j]-vmin[j];
  219.         If((stdv<>0.0) and (nc[j]>2.0)) Then 
  220.            skw:=nc[j]*mm3[j]/((nc[j]-1.0)*(nc[j]-2.0)*vrnc*stdv)
  221.            Else skw:=0.0;
  222.         If((vrnc<>0.0) And (nc[j]>3.0)) Then krt:=(nc[j]*(nc[j]+1.0)*
  223.            mm4[j]-Sqr(mm2[j])*(nc[j]-1)*3.0)/((nc[j]-1.0)*(nc[j]-2.0)*
  224.            (nc[j]-3.0)*Sqr(vrnc))
  225.            Else krt:=0.0;
  226.         Writeln(ofile,'Mean:    ',mu:13:5,'   Std. Error:',ster:13:5);
  227.         Writeln(ofile,'Variance:',vrnc:13:5,'   Std. Dev:  ',stdv:13:5);
  228.         Writeln(ofile,'Skewness:',skw:13:5,'   Kurtosis:  ',krt:13:5);
  229.         Writeln(ofile,'Max:     ',vmax[j]:13:5,'   Min:       ',vmin[j]:13:5);
  230.         Writeln(ofile,'Range:   ',rng:13:5,'   Cases:     ',nc[j]:13:5);
  231.         Writeln(ofile);
  232.         (* construct frequency histogram *)
  233.         If(gp<>0) Then
  234.          Begin
  235.          Write('- Press any key to continue -');
  236.          While Not KeyPressed Do;
  237.          End;
  238.         If(gp=0) Then
  239.          Begin
  240.          Writeln('Histograms take another pass at the data, do you want');
  241.          Write('one for variable ',varn[j],'? '); Readln(yn);
  242.          If((Copy(yn,1,1)='Y') Or (Copy(yn,1,1)='y')) Then
  243.           Begin
  244.             rng:=ot*20.0; If rng>nc[j] Then rng:=nc[j];
  245.             Write('Intervals in histogram (recommend',Trunc(rng):3,')? ');
  246.             Readln(rng); Writeln;
  247.             histogram(vmax[j],vmin[j],rng,sel[j],varn[j],miss[j],
  248.                       dfile,ofile);
  249.             wait(ot);
  250.             End;
  251.           End;
  252.         End;
  253.       ClrScr;
  254.       End;
  255.    If Not EOF(dfile) Then Begin fg:=1; Goto initer; End;
  256.    Close(dfile); Close(ofile);
  257.    Assign(dfile,'MAPSTAT.COM'); Execute(dfile);
  258. End.
  259.  ',rng:13:5,'   Cases:     ',n