% LMcs110206_4th.m
% Fitting a logistic regression
% Written 11/10 by Eugene.Gallagher@umb.edu
% Revised 2/20/2011. 
DATA=[480 0.3;690 4.6;900 15.6;1100 33.4;1320 44.4;1530 45.7];
X=DATA(:,1);y=DATA(:,2)/100;
[b,dev,stats]=glmfit(X,y,'binomial')
% [yhat,dylo,dyhi] = glmval(b,X,'logit',stats)
yhat = glmval(b,X,'logit')
plot(X,y,'o',X,yhat,'-b','LineWidth',2);
xlabel('SAT Score');ylabel('Graduation Rate (%)')
figure(gcf);pause

% The fit looks really awful
% Regress ln (L-y/y) vs X as described in Larsen & Marx (2006) p 669
Y=log((48-DATA(:,2))./DATA(:,2));
[r,c]=size(Y);
X=[ones(r,1) DATA(:,1)];
B=X\Y;
fprintf('\nThe Y intercept = %8.3f and slope = %7.5f\n',B(1),B(2));
Yest=X*B;
Resid=Y-Yest;
% sort by the independent variable so the lines are smooth
[Xs,k]=sort(X(:,2));
% sort the estimated values using the index vector k:
Ys=Yest(k);
plot(X(:,2),Y,'o',Xs,Ys,'m','LineWidth',2,...
                       'MarkerEdgeColor','k',...
                       'MarkerFaceColor','b',...
                       'MarkerSize',8) 
xlabel('SAT','FontSize',20),ylabel('Graduation %','FontSize',20)
title('Figure 11.2.10a','FontSize',22);
figure(gcf)
pause
% Generated a fine-scale estimated curve
xintrp=[1:0.2:1600]';
[r,c]=size(xintrp);
xintrp=[ones(r,1) xintrp];
yestintrp=xintrp*B;
plot(DATA(:,1),DATA(:,2),'o',...
    xintrp(:,2),48./ (1+exp(B(1)+B(2).*xintrp(:,2) )),'--r',...
                        'LineWidth',2,...
                       'MarkerEdgeColor','k',...
                       'MarkerFaceColor','b',...
                       'MarkerSize',8) 
   xlabel('SAT','FontSize',20),ylabel('Graduation %','FontSize',20)
   title('Figure 11.2.10','FontSize',22)
axis([0 1650 0 55])
figure(gcf)
pause
plot(Yest,Resid,'ok');
xlabel('Expected value'),ylabel('Residual');
grid;figure(gcf)
