Ellie regression

############################
###GUM2011AllPositivity#####
############################
summary(dat$totalPositivity2011GUM)
hist(dat$totalPositivity2011GUM)

modpos<-glm(totalPositivity2011GUM~ratioWhite2011+ratioBlack2011+ratioAsian2011+
	conception.rate.per.1000women.Under18+MarriedCivilCohabitCouples2011+SingleParents2011+
	ratioNoQual_2011+ratioMales25to34_2011+RatioFemales16to24_2011,
	family=binomial,data=dat,weights=sum.NtestsGUM2011All)
summary(modpos)
theta<-modpos$deviance/modpos$df.residual
theta
modpos<-glm(totalPositivity2011GUM~conception.rate.per.1000women.Under18+MarriedCivilCohabitCouples2011+
	ratioNoQual_2011,
	family=quasibinomial,data=dat,weights=sum.NtestsGUM2011All)
drop1(modpos,test="F")
summary.lm(modpos)

#modpos<-glm(totalPositivity2011GUM~ratioBlack2011+
#	conception.rate.per.1000women.Under18+MarriedCivilCohabitCouples2011+SingleParents2011+
#	ratioNoQual_2011,
#	family=binomial,data=dat,weights=sum.NtestsGUM2011All)
#drop1(modpos,test="Chi")

summary(modpos)
par(mfrow=c(2,4))
plot(modpos)
hist(modpos$resid)
plot(modpos$resid~dat$ratioBlack2011)
plot(modpos$resid~dat$conception.rate.per.1000women.Under18)
plot(modpos$resid~dat$MarriedCivilCohabitCouples2011)
plot(modpos$resid~dat$SingleParents2011)
plot(modpos$resid~dat$ratioNoQual_2011)

summary(dat$ratioNoQual_2011)
(((
par(mfrow=c(2,3))
MyDat<-data.frame(ratioBlack2011=seq(0.000419,0.271600,0.0008253317),
	conception.rate.per.1000women.Under18=mean(dat$conception.rate.per.1000women.Under18),
	MarriedCivilCohabitCouples2011=mean(dat$MarriedCivilCohabitCouples2011),
	SingleParents2011=mean(dat$SingleParents2011),
	ratioNoQual_2011=mean(dat$ratioNoQual_2011))
head(MyDat)
P1<-predict(modpos,newdata=MyDat,type="link",se=TRUE)
plot(MyDat$ratioBlack2011, exp(P1$fit)/(1+exp(P1$fit)),
	type="l",
	xlab="Ratio of Black individuals",ylab="Probability of chlamydia infection") lines(MyDat$ratioBlack2011, exp(P1$fit+1.96*P1$se.fit)/
	(1+exp(P1$fit+1.96*P1$se.fit)),lty=2)
lines(MyDat$ratioBlack2011, exp(P1$fit-1.96*P1$se.fit)/
	(1+exp(P1$fit-1.96*P1$se.fit)),lty=2)
points(dat$ratioBlack2011,dat$totalPositivity2011GUM)

MyDat2<-data.frame(ratioBlack2011=mean(dat$ratioBlack2011),
	conception.rate.per.1000women.Under18=seq(9,58.1,0.1515418),
	MarriedCivilCohabitCouples2011=mean(dat$MarriedCivilCohabitCouples2011),
	SingleParents2011=mean(dat$SingleParents2011),
	ratioNoQual_2011=mean(dat$ratioNoQual_2011))
P2<-predict(modpos,newdata=MyDat2,type="link",se=TRUE)
plot(MyDat2$conception.rate.per.1000women.Under18, exp(P2$fit)/(1+exp(P2$fit)),
	type="l",
	xlab="Conception Rate U18 per1000women",ylab="Probability of chlamydia infection",
	main="GUMCAD positivity data for all patients") lines(MyDat2$conception.rate.per.1000women.Under18, exp(P2$fit+1.96*P2$se.fit)/
	(1+exp(P2$fit+1.96*P2$se.fit)),lty=2)
lines(MyDat2$conception.rate.per.1000women.Under18, exp(P2$fit-1.96*P2$se.fit)/
	(1+exp(P2$fit-1.96*P2$se.fit)),lty=2)
points(dat$conception.rate.per.1000women.Under18,dat$totalPositivity2011GUM)

MyDat3<-data.frame(ratioBlack2011=mean(dat$ratioBlack2011),
	conception.rate.per.1000women.Under18=mean(dat$conception.rate.per.1000women.Under18),
	MarriedCivilCohabitCouples2011=seq(0.4423,0.8007,0.001106162),
	SingleParents2011=mean(dat$SingleParents2011),
	ratioNoQual_2011=mean(dat$ratioNoQual_2011))
P3<-predict(modpos,newdata=MyDat3,type="link",se=TRUE)
plot(MyDat3$MarriedCivilCohabitCouples2011, exp(P3$fit)/(1+exp(P3$fit)),
	type="l",
	xlab="Ratio of individuals in live-in stable partnerships",ylab="Probability of chlamydia infection") lines(MyDat3$MarriedCivilCohabitCouples2011, exp(P3$fit+1.96*P3$se.fit)/
	(1+exp(P3$fit+1.96*P3$se.fit)),lty=2)
lines(MyDat3$MarriedCivilCohabitCouples2011, exp(P3$fit-1.96*P3$se.fit)/
	(1+exp(P3$fit-1.96*P3$se.fit)),lty=2)
points(dat$MarriedCivilCohabitCouples2011,dat$totalPositivity2011GUM)

MyDat4<-data.frame(ratioBlack2011=mean(dat$ratioBlack2011),
	conception.rate.per.1000women.Under18=mean(dat$conception.rate.per.1000women.Under18),
	MarriedCivilCohabitCouples2011=mean(dat$MarriedCivilCohabitCouples2011),
	SingleParents2011=seq(0.04828,0.19910,0.0004654894),
	ratioNoQual_2011=mean(dat$ratioNoQual_2011))
P4<-predict(modpos,newdata=MyDat4,type="link",se=TRUE)
plot(MyDat4$SingleParents2011, exp(P4$fit)/(1+exp(P4$fit)),
	type="l",
	xlab="Ratio of single parents",ylab="Probability of chlamydia infection") lines(MyDat4$SingleParents2011, exp(P4$fit+1.96*P4$se.fit)/
	(1+exp(P4$fit+1.96*P4$se.fit)),lty=2)
lines(MyDat4$SingleParents2011, exp(P4$fit-1.96*P4$se.fit)/
	(1+exp(P4$fit-1.96*P4$se.fit)),lty=2)
points(dat$SingleParents2011,dat$totalPositivity2011GUM)

MyDat5<-data.frame(ratioBlack2011=mean(dat$ratioBlack2011),
	conception.rate.per.1000women.Under18=mean(dat$conception.rate.per.1000women.Under18),
	MarriedCivilCohabitCouples2011=mean(dat$MarriedCivilCohabitCouples2011),
	SingleParents2011=mean(dat$SingleParents2011),
	ratioNoQual_2011=seq(0.06721,0.3517,0.0008780472))
P5<-predict(modpos,newdata=MyDat5,type="link",se=TRUE)
plot(MyDat5$ratioNoQual_2011, exp(P5$fit)/(1+exp(P5$fit)),
	type="l",
	xlab="Ratio of people without qualifications",ylab="Probability of chlamydia infection") lines(MyDat5$ratioNoQual_2011, exp(P5$fit+1.96*P5$se.fit)/
	(1+exp(P5$fit+1.96*P5$se.fit)),lty=2)
lines(MyDat5$ratioNoQual_2011, exp(P5$fit-1.96*P5$se.fit)/
	(1+exp(P5$fit-1.96*P5$se.fit)),lty=2)
points(dat$ratioNoQual_2011,dat$totalPositivity2011GUM)