%include catdata(titanic);

*-- Handle the many zero cells.  Note that GENMOD treats all such cells
    as structural zeros.  This causes some of the models to fail to
	converge, so we treat them as sampling zeros by adding a small positive
	number to each frequency;
*-- (We should probably treat Children among the Crew as structural zeros);
data titanic;
	set titanic;
	count = count+.5;

options ls=90;

*-- Because Survival is considered the response, we must
    include the high-order term among the explanatory variables;
title '[AGC][S]: Baseline model'; 
proc genmod data=titanic;
	class age sex class survive;
	model count=age|sex|class survive  / dist=poisson type3 obstats;
	ods output obstats=obstats;
%mosaic(data=obstats, vorder=class sex age survive,
	resid=reschi, title=[AGC][S]: Baseline model)

*-- 'Main effects model for survival';
title '[AGC][AS][GS][CS]: Associations of age, sex, class with survive';
proc genmod data=titanic;
	class age sex class survive;
	model count=age|sex|class 
				age*survive sex*survive class*survive
				 / dist=poisson type3 obstats;
	ods output obstats=obstats;
%mosaic(data=obstats, vorder=class sex age survive,
	resid=reschi, title=[AGC][AS][GS][CS])

title '[AGC][CS][GAS]';
proc genmod data=titanic;
	class age sex class survive;
	model count=age|sex|class 
				age|sex|survive class*survive
				 / dist=poisson type3 obstats;
	ods output obstats=obstats;
%mosaic(data=obstats, vorder=class sex age survive,
	resid=reschi, title=[AGC][CS][GAS])


title '[AGC]CGS][GAS]';
proc genmod data=titanic;
	class age sex class survive;
	model count=age|sex|class 
				age|sex|survive class|sex|survive
				 / dist=poisson type3 obstats;
	ods output obstats=obstats;
%mosaic(data=obstats, vorder=class sex age survive,
	resid=reschi, title=[AGC][CGS][GAS])


%genmodsum(data=titanic, yvariable=count,
	class=age sex class survive,
	model1 =age|sex|class survive / [ASC][S],
	model2 =age|sex|class 
				age*survive sex*survive class*survive/ [AGC][AS][GS][CS],
	model3 = age|sex|class 
				age|sex|survive class*survive/ [AGC]CS][GAS],
	model4 = age|sex|class 
				age|sex|survive class|sex|survive/ [AGC]CGS][GAS],
	distribution=poisson)  
;