%MATLAB SCRIPT FILE s322eg7.m  */
%
%    EXAMPLE 7 FOR STATISTICS 322
%    POSTED 9/17/96
%    EXPLORES Q-Q AND P-P PLOTS

ipar = 1 ;     %  1,2

autop = 1 ;    %  1 TO SAVE PLOTS AS POSTSCRIPT FILES


disp('Running MATLAB script s322eg6.m') ;


if ipar == 1 ;    %  THEN DO PART 1: Q-Q PLOTS

    npgrid = 20 ;
    pdel = 1 / (npgrid + 1) ;
  pgrid = linspace(pdel,1-pdel,npgrid) ;
          %  EQUALLY SPACED GRID IN SENSE OF MEAN OF UNIF. ORDER STATS.

  for ieg = 1:4 ;

    if ieg == 1 ;
      subplot(2,2,1) ;
        qx = norminv(pgrid,0,1) ;
          %  VECTOR OF N(0,1) QUANTILES
        qy = norminv(pgrid,5,2) ;
          %  VECTOR OF N(5,4) QUANTILES
          %  CAREFUL:  THAT LAST ARGUMENT IS S.D. NOT VAR!
      plot(qx,qy,'+',qx,qy,'-') ;
          %  USE TWO VERSIONS TO GET PLUS SIGNS CONNECTED BY LINES
        title('Q-Q Plot of N(0,1) vs. N(5,4)') ;
        xlabel('N(0,1)') ;
        ylabel('N(5,4)') ;
        axis('square') ;
        axis('equal') ;

    elseif ieg == 2 ;
      subplot(2,2,2) ;
      qx = norminv(pgrid,0,1) ;
      qy = norminv(pgrid,-1,2/3) ;
      plot(qx,qy,'+',qx,qy,'-') ;
        title('Q-Q Plot of N(0,1) vs. N(-1,4/9)') ;
        xlabel('N(0,1)') ;
        ylabel('N(-1,4/9)') ;
        axis('square') ;
          %  THIS CONTRAINS THE PLOT TO LIE IN A SQUARE REGION,
          %  VS. THE USUAL "FILLING UP AS MUCH SPACE AS POSSIBLE"
        axis('equal') ;
          %  THIS MAKES LENGTHS ON THE AXES BE IN THE SAME UNITS

          %  THE ABOVE TWO COMMANDS ARE USEFUL WHEN "VERTICAL LENGTH"
          %  SHOULD HAVE THE SAME VISUAL CONTENT AS "HORIZONTAL LENTH"

    elseif ieg == 3 ;
      subplot(2,2,3) ;
      qy = tan(pi * (pgrid - 1/2)) ;
          %  RECALL CAUCHY CDF IS 1/2 + ARCTAN(X)/PI
        q3c = tan(pi * (.75 - 1/2)) ;
          %  THIRD QUARTILE OF CAUCHY
        q3n = norminv(.75,0,1) ;
          %  THIRD QUARTILE OF N(0,1)
      qx = norminv(pgrid,0,q3c/q3n) ;
      plot(qx,qy,'+',qx,qy,'-') ;
        title('Q-Q Plot of N(0,q3c/qcn) vs. Cauchy') ;
        xlabel('N(0,q3c/q3n)') ;
        ylabel('Cauchy') ;
        axis('square') ;
        axis('equal') ;

    elseif ieg == 4 ;
      subplot(2,2,4) ;
      qx = norminv(pgrid,1/2,sqrt(1/12)) ;
      qy = pgrid ;
          %  UNIFORM(0,1) C.D.F. 
      plot(qx,qy,'+',qx,qy,'-') ;
        title('Q-Q Plot of N(1/2,1/12) vs. Uniform') ;
        xlabel('N(1/2,1/12)') ;
        ylabel('Uniform') ;
        axis('square') ;
        axis('equal') ;

    end ;

  end ;

  if autop == 1 ;
    print s322eg7a.ps -dps ;
  end ;
          %  PUTTING THIS INSIDE THIS if BLOCK ALLOWS ONE TO CHOOSE
          %  WHETHER OR NOT TO SAVE TO FILE (USEFUL WHILE DEBUGGING).


elseif ipar == 2 ;    %  THEN DO PART 2: P-P PLOTS


end ;

