Regression Models with Spatially Misaligned Data

MISSING IMAGE

Material Information

Title:
Regression Models with Spatially Misaligned Data
Physical Description:
1 online resource (148 p.)
Language:
english
Creator:
Lopiano, Kenneth K
Publisher:
University of Florida
Place of Publication:
Gainesville, Fla.
Publication Date:

Thesis/Dissertation Information

Degree:
Doctorate ( Ph.D.)
Degree Grantor:
University of Florida
Degree Disciplines:
Statistics
Committee Chair:
Young, Linda
Committee Members:
Daniels, Michael J
Casella, George
Mcintyre, Lauren

Subjects

Subjects / Keywords:
generalized -- kriging -- linear -- model -- regression -- spatial
Statistics -- Dissertations, Academic -- UF
Genre:
Statistics thesis, Ph.D.
bibliography   ( marcgt )
theses   ( marcgt )
government publication (state, provincial, terriorial, dependent)   ( marcgt )
born-digital   ( sobekcm )
Electronic Thesis or Dissertation

Notes

Abstract:
In many applications, the response Y and the covariate X are associated with geographical units (e.g., a point in space or an areal unit). In practice, the geographical units associated with the observed response values and observed covariate values may be different. As a result, the data are said to be spatially misaligned. For example, if Y is observed for areal units and X is observed at points in space, a point-to-areal misalignment problem arises. As another example, if Y is observed at one set of points in space and X is observed at a different set of points in space, a point-to-point misalignment problem arises.  Under certain modeling assumptions, the conditional distribution of the covariate at the geographical units associated with the response given the observed values of the covariate can be derived. The mean of this distribution is referred to as the kriging mean. The kriging mean has a Berkson error relationship with the true covariate. Naive inferences that use the kriging mean in subsequent regression models and ignore the additional uncertainty (i.e., the variability associated with conditional distribution) can result in biased effect estimates and understating the uncertainty associated with these estimates.  We consider the Berkson error problem in generalized linear models and develop methods to estimate regression parameters and associated measures of uncertainty. Using the spatial misalignment problem as a motivating example throughout, we first review and compare existing methods developed to address the spatial misalignment problem in normal linear models. Then we show how methods used to conduct inference in linear mixed models and generalized linear mixed models can be adapted to conduct inference in generalized linear models with Berkson error. The properties of the inference procedures are assessed theoretically and via simulation and illustrated using examples of spatial misalignment.
General Note:
In the series University of Florida Digital Collections.
General Note:
Includes vita.
Bibliography:
Includes bibliographical references.
Source of Description:
Description based on online resource; title from PDF title page.
Source of Description:
This bibliographic record is available under the Creative Commons CC0 public domain dedication. The University of Florida Libraries, as creator of this bibliographic record, has waived all rights to it worldwide under copyright law, including all related and neighboring rights, to the extent allowed by law.
Statement of Responsibility:
by Kenneth K Lopiano.
Thesis:
Thesis (Ph.D.)--University of Florida, 2012.
Local:
Adviser: Young, Linda.
Electronic Access:
RESTRICTED TO UF STUDENTS, STAFF, FACULTY, AND ON-CAMPUS USE UNTIL 2014-08-31

Record Information

Source Institution:
UFRGP
Rights Management:
Applicable rights reserved.
Classification:
lcc - LD1780 2012
System ID:
UFE0044472:00001


This item is only available as the following downloads:


Full Text

PAGE 1

REGRESSIONMODELSWITHSPATIALLYMISALIGNEDDATAByKENNETHKYLELOPIANOADISSERTATIONPRESENTEDTOTHEGRADUATESCHOOLOFTHEUNIVERSITYOFFLORIDAINPARTIALFULFILLMENTOFTHEREQUIREMENTSFORTHEDEGREEOFDOCTOROFPHILOSOPHYUNIVERSITYOFFLORIDA2012

PAGE 2

c2012KennethKyleLopiano 2

PAGE 3

ToDr.CasellaYourunwaveringpassionhasinspiredsomanyIfyouhaveanappleandIhaveanappleandweexchangetheseapplesthenyouandIwillstilleachhaveoneapple.ButifyouhaveanideaandIhaveanideaandweexchangetheseideas,theneachofuswillhavetwoideas.-GeorgeBernardShaw 3

PAGE 4

ACKNOWLEDGMENTS Firstandforemost,Iwouldlikeacknowledgemyadvisor,Dr.LindaYoung.Dr.Young,youhavetaughtmesomuchaboutbeingagreatstatisticianandagreatperson.IwillalwaysrememberoneofourrstexperiencestogetheratNISSduringabrainstormingsession.ItwaswhileworkingwithyouatNISSthatIfellinlovewithresearch.Sincethen,youhaveinvestedsomuchofyourtimeandenergytoensurethatIwouldbeabletogrowtomyfullpotential.Thankyousomuchforyourpatienceandguidanceduringthepastthreeyears.Thankyouforteachingmehowtoasktherightquestionsandteachingmethatresearchistwostepsforwardandonestepbackward.AsIcompletemydissertationunderyourguidance,IknowthatthisisoneofmanystepsforwardIwilltakebecauseofyou.Iamtrulyhonoredtobeyourstudent.Second,IwouldliketoacknowledgealloftheprofessorsattheUniversityofFloridathatprovidedmewiththeskillsneededtobewhereIamtoday.Inparticular,IwouldliketothankDr.Presnell,Dr.Doss,andDr.Casellaforpushingmetobethebeststudentpossibleduringmysecondyear,Dr.McIntyreforsomanyinvaluablelessons,andDr.Gotwayforhelpingmeinsomanyways.Dr.McIntyre,thankyouforallthelessonsaboutbeingagoodscientistandgivingmetheopportunitytocollaboratewithyou.Dr.Gotway,Ihavehadsomuchfunworkingwithyouoverthepastfewyears.OurconversationsaboutstatisticsandlifearepricelessandIamsogratefultohavehadtheopportunitytolearnsomuchfromyou.Dr.Presnell,thankyouforthecountlesshoursyouspentworkingwithmeonprobabilitytheoryandsincerelytakinganinterestinmyfuture.Iwouldnotbeherewithoutyourguidanceandmotivation.Dr.Doss,thankyouforyourhonestyandencouragement.Dr.Casella,itwasanhonortolearnsomuchfromyouaboutstatisticsandlife.HearingyourlaughterechothroughthehallsofMcCartyalwaysremindedmetoenjoytheexperienceandhavefun.Seeingthephotographsofyourchildrenonyourcomputerandhearinghowproudyouwereofthemshowedmehowimportantitistobalance 4

PAGE 5

familyandwork.YourpassionforstatisticsandlifeweretrulyinspiringandIamsohappytohavehadtheopportunitytolearnfromyouinsomanyways.Finally,Iwouldliketoacknowledgemyfamilyandfriends.Tomybrother,Mickey,thankyousomuchforbeingtherstpersontoteachmemathematics.TomyDadandAimee,thankyouforalwaysbeingtheretosupportmeinsomanyways.Tomygirlfriend,Jeannie,thankyouforalwayssupportingmeandencouragingme.Yoursupportduringthenalstretchhelpedmeinmorewaysthanyouknow.Finally,tomyfriendsNateandMike.Mike,thankyouforalwaysbeingtheretolistenandprovidingmewiththecalmingadvicethatIneeded.Nate,wetraveledthisroadtogether.IhadsomuchfunworkingwithyouonprobabilitytheoryproblemsinGrifn-Floyd,studyingtogetherforthesecondyearexam,andsittinginmylivingroomdiscussingproblemsandideas.You'reagreatfriend.Icannotthankyouenoughforthelessonsyou'vetaughtmeoverthepastfewyears. 5

PAGE 6

TABLEOFCONTENTS page ACKNOWLEDGMENTS .................................. 4 LISTOFTABLES ...................................... 9 LISTOFFIGURES ..................................... 10 ABSTRACT ......................................... 11 CHAPTER 1INTRODUCTION ................................... 13 2ACOMPARISONOFERRORSINVARIABLESMETHODSFORUSEINREGRESSIONMODELSWITHSPATIALLYMISALIGNEDDATA ................. 16 2.1Background ................................... 16 2.2StatisticalModelsAndMethods ........................ 19 2.2.1TraditionalKrigeAndRegress(KR) .................. 22 2.2.2TraditionalKrigeAndRegressWithGeneralCovarianceStructure(KRGC) ................................. 22 2.2.3AdjustedKrigeAndRegressApproaches .............. 22 2.2.3.1Naiveadjustedkrigeandregress(NAKR) ......... 23 2.2.3.2Variance-correctedadjustedKR(VCAKR) ......... 23 2.2.4ParametricBootstrap(PB) ....................... 24 2.2.5ParameterBootstrap(PRB) ...................... 26 2.2.6PartialParametricBootstrap(PPB) .................. 26 2.3TheStudyAndSupportingData ....................... 27 2.4SimulationStudyForMisalignedPoints ................... 28 2.4.1Scenario1 ................................ 29 2.4.2Scenario2 ................................ 31 2.5MovingFromPointstoAreas ......................... 32 2.5.1ResultsForTheFloridaData ..................... 33 2.5.2FloridaGridSimulation ......................... 34 2.5.3Scenario2Simulation ......................... 34 2.6ConcludingRemarks .............................. 35 2.7Support ..................................... 38 3SPATIALCOVARIANCEESTIMATIONINSPATIALLY-MISALIGNEDREGRESSIONMODELSWITHBERKSONERROR ........................ 42 3.1Background ................................... 42 3.2ModelingFramework .............................. 46 3.3ProposedMethodology ............................ 47 3.3.1GLSEstimation ............................. 49 3.3.1.1Method-of-moments ..................... 52 6

PAGE 7

3.3.1.2Restrictedmaximumlikelihoodbasedmethods ...... 53 3.3.2EstimatingTheBerksonParameters ................. 55 3.4Applications ................................... 57 3.4.1EPAEnvironmentalMonitoringandAssessmentProgram ..... 57 3.4.2CDCEnvironmentalPublicHealthTrackingProgram ........ 59 3.5Simulation .................................... 62 3.5.1Point-to-PointMisalignment ...................... 62 3.5.2Point-to-ArealMisalignment ...................... 63 3.6ConcludingRemarks .............................. 65 4BERKSONERRORINGENERALIZEDLINEARMODELS ........... 75 4.1BackgroundInformation ............................ 75 4.2ProposedMethodology ............................ 80 4.3ComparingNaiveandAdjustedEstimators ................. 84 4.4SimulationStudy ................................ 89 4.4.1StudyResults .............................. 90 4.4.1.1Poisson ............................ 90 4.4.1.2Binomial ........................... 91 4.4.1.3Normal ............................ 92 4.5Application ................................... 93 4.5.1Point-to-ArealMisalignment ...................... 93 4.5.2Point-to-PointMisalignment ...................... 96 4.6ConcludingRemarks .............................. 98 5ASYMPTOTICPROPERTIESOFTHEMOMESTIMATOR ........... 110 6ASYMPTOTICPROPERTIESOFPSEUDO-LIKELIHOODESTIMATORSINLINEARREGRESSIONMODELS ......................... 115 6.1ModelFormulation ............................... 115 6.2GLS ....................................... 115 6.3EGLS ...................................... 116 6.4OLS ....................................... 116 6.5REML ...................................... 118 6.6PseudoRestrictedLikelihood ......................... 120 6.7EGLSEstimator ................................ 123 7FUTUREWORK ................................... 124 7.1Spatio-TemporalMisalignment ........................ 124 7.2ClassicalMeasurementError ......................... 124 7.2.1PoissonRegressionwithCME ..................... 124 7.2.2Simulation ................................ 125 7.3BerksonError .................................. 125 7.3.1PoissonRegression:UnknownBerksonParameter ......... 126 7.3.2Simulation ................................ 126 7

PAGE 8

7.4TablesandFigures ............................... 127 APPENDIX:SUPPLEMENTARYMATERIALS ...................... 128 A.1ShowingQ2=Q3 ................................ 128 A.2EigenvaluesofInducedCovarianceMatrices ................ 128 A.3ExpectedHealthCounts ............................ 129 REFERENCES ....................................... 143 BIOGRAPHICALSKETCH ................................ 147 8

PAGE 9

LISTOFTABLES Table page 2-1Scenario1SimulationResults .......................... 39 2-2Scenario2SimulationResults .......................... 39 2-3Estimatesof^1andtheStandardErrorof^1forFlorida ............ 39 2-4Scenario2SimulationResults .......................... 40 2-5Average^1andVariance^2^1oftheEstimatesof1andAverageoftheEstimatedVariances2^1of^1fortheSimulatedDatasetsintheFloridaGridSimulation 40 3-1EMAPApplication:Estimatedregressioncoefcientandassociatedstandarderrorsbasedondifferentmethods. ........................ 68 3-2EPHTApplication:Estimatedregressioncoefcientsandmeasuresofuncertaintybasedondifferentmethods ............................ 68 3-3EPHTApplication:Estimatesofregressionerrorbasedondifferentmethods 68 3-4Point-to-PointSimulationResults ......................... 69 3-5Point-to-ArealSimulationResults ......................... 69 4-1PoissonResponseSimulationResults ...................... 102 4-2BinomialResponseSimulationResults ...................... 103 4-3NaiveInference ................................... 104 4-4NormalResponseSimulationResults ...................... 104 7-1ClassicalMeasurementErrorSimulationResults ................ 127 7-2BerksonErrorSimulationResults ......................... 127 9

PAGE 10

LISTOFFIGURES Figure page 2-1CoarsegridlaidoverthestateofFlorida ...................... 40 2-2AirmonitorlocationsinthestateofFlorida ..................... 41 3-1ObservationLocations ................................ 70 3-2Fittedregressionlinesfordifferentmethods. .................... 71 3-3MonitorLocationsandCountiesintheStateofFlorida .............. 72 3-4LocationsofPointsforSimulationStudy ...................... 73 3-5GridLocationsandCountiesintheStateofFlorida ................ 74 4-1PoissonResponseNaiveCoverage ........................ 105 4-2PoissonResponsePQL-basedCoverage ..................... 106 4-3PoissonResponseNaive:PQL-basedratio ..................... 107 4-4BinomialResponseNaive:PQL-basedratio .................... 108 4-5PM2.5monitorsinFlorida .............................. 109 A-1Simulationlocations ................................. 131 A-2Correlationfunctions ................................. 132 A-3Poissonmeans .................................... 133 A-4Binomialmeans ................................... 134 A-5PoissonResponseNaiveVariance ......................... 135 A-6PoissonResponsePQL-basedVariance ...................... 136 A-7BinomialResponseNaiveVariance ......................... 137 A-8BinomialResponseNaiveCoverage ........................ 138 A-9BinomialResponsePQL-basedVariance ..................... 139 A-10BinomialResponsePQL-basedCoverage ..................... 140 A-11NASAsatelliteimagefromApril29th,2007 .................... 141 A-12NASAsatelliteimagefromApril30th,2007 .................... 141 A-13FinegridlaidoverFlorida .............................. 142 10

PAGE 11

AbstractofDissertationPresentedtotheGraduateSchooloftheUniversityofFloridainPartialFulllmentoftheRequirementsfortheDegreeofDoctorofPhilosophyREGRESSIONMODELSWITHSPATIALLYMISALIGNEDDATAByKennethKyleLopianoAugust2012Chair:Dr.LindaJ.YoungMajor:StatisticsInmanyapplications,theresponseYandthecovariateXareassociatedwithgeographicalunits(e.g.,apointinspaceoranarealunit).Inpractice,thegeographicalunitsassociatedwiththeobservedresponsevaluesandobservedcovariatevaluesmaybedifferent.Asaresult,thedataaresaidtobespatiallymisaligned.Forexample,ifYisobservedforarealunitsandXisobservedatpointsinspace,apoint-to-arealmisalignmentproblemarises.Asanotherexample,ifYisobservedatonesetofpointsinspaceandXisobservedatadifferentsetofpointsinspace,apoint-to-pointmisalignmentproblemarises.Undercertainmodelingassumptions,theconditionaldistributionofthecovariateatthegeographicalunitsassociatedwiththeresponsegiventheobservedvaluesofthecovariatecanbederived.Themeanofthisdistributionisreferredtoasthekrigingmean.ThekrigingmeanhasaBerksonerrorrelationshipwiththetruecovariate.Naiveinferencesthatusethekrigingmeaninsubsequentregressionmodelsandignoretheadditionaluncertainty(i.e.,thevariabilityassociatedwithconditionaldistribution)canresultinbiasedeffectestimatesandunderstatingtheuncertaintyassociatedwiththeseestimates.WeconsidertheBerksonerrorproblemingeneralizedlinearmodelsanddevelopmethodstoestimateregressionparametersandassociatedmeasuresofuncertainty.Usingthespatialmisalignmentproblemasamotivatingexamplethroughout,werst 11

PAGE 12

reviewandcompareexistingmethodsdevelopedtoaddressthespatialmisalignmentprobleminnormallinearmodels.ThenweshowhowmethodsusedtoconductinferenceinlinearmixedmodelsandgeneralizedlinearmixedmodelscanbeadaptedtoconductinferenceingeneralizedlinearmodelswithBerksonerror.Thepropertiesoftheinferenceproceduresareassessedtheoreticallyandviasimulationandillustratedusingexamplesofspatialmisalignment. 12

PAGE 13

CHAPTER1INTRODUCTIONResearchersroutinelyusegeneralizedlinearmodelstoassesstherelationshipbetweenaresponseYandacovariateX.Morespecically,researchersareofteninterestedinassessingtherelationshipsamongvariableswithexplicitlyorimplicitlydenedgeographicalunits.Forexample,inenvironmentalpublichealthstudies,therelationshipbetweenthebirthweightofaninfant,linkedgeographicallytotheresidenceofitsmother,andtheconcentrationofozoneorparticulatematter(PM),linkedgeographicallytoamonitorlocation,mightbeofinterest.Regardlessoftheapplication,thevariablesareoftenobservedatdifferentlocationsoraggregatedoverdifferentgeographicalunits.Suchdataaresaidtobespatiallymisaligned.Forconcreteness,supposetheresponsevariableYisobservedatnspatiallocationssu,andtheexplanatoryvariableXisobservedatnodifferentspatiallocationsso.SupposetherelationshipbetweenY(su)andX(su)followsageneralizedlinearmodelY(su)ExponentialFamily(mean=,scale=)g()=01n1+1X(su), (1)wheregisamonotoniclinkfunction.Beforeinferencecanbeconducted,thedisparatedatasetsmustbealigned.Inpractice,ajointdistributionisassumedfortheunobservedandobservedcovariatevalues.LetX(su)betheunobservedXvaluesatthelocationssuandX(so)betheobservedXvaluesatthelocationsso.SupposethejointdistributionofX(su)andX(so)is 264X(su)X(so)375N0B@264C(su)C(so)375,X(x)=264X,uuX,uoX,ouX,oo3751CA,(1)whereC(su)andC(so)arenkandnokmatricesofobservedcovariatesatthelocationsinsuandsorespectively,isanunknownvectorkparameters,xisan 13

PAGE 14

unknownvectoroflparametersthatdenethecovariancematrixXwhichispartitionedintothecovariancematricesX,uu,thecovarianceamongtheunobservedlocations,X,uo,thecovarianceamongtheunobservedandobservedlocations,andX,oo,thecovarianceamongtheobservedlocations.Forexample,xmaybeparametersthatdeneanexponentialcovariancestructure.Theconditional(kriging)distributionofX(su)jX(so)is X(su)jX(so),,xN(W=C(su)+X,uo)]TJ /F5 7.97 Tf 6.59 0 Td[(1X,oo(X(so))]TJ /F10 11.955 Tf 11.96 0 Td[(C(so))x=X,uu)]TJ /F10 11.955 Tf 11.96 0 Td[(X,uo)]TJ /F5 7.97 Tf 6.59 0 Td[(1X,ooX,ou). (1) NoticethekrigingmeanhasaBerksonerrorrelationshipwiththeunobservedcovariate, X(su)jX(so),,x=W+x,(1)wherexN(0,x).Generallyspeaking,aBerksonerrorrelationshipisameasurementerrormodelthattakestheform X=W+x.(1)wherexisameanzerorandomvariableindependentofWandY.BecauseofthisBerksonerrorrelationship,Y(su)ExponentialFamily(mean=,scale=)g()=01n1+1W+1xxN(0,x). (1)TheBerksonerrorrelationshipinthecontextofspatialmisalignmenthasbeenidentiedanddiscussedbyseveralresearchers( Gryparisetal. ( 2009 ); Lopianoetal. ( 2011 ); Madsenetal. ( 2008 ); Pacioreketal. ( 2009 ); Szpiroetal. ( 2011 )).BecauseBerksonerrorcancauseestimatesof1anditsassociatedmeasuresofuncertaintytobebiased,itisofinteresttodevelopinferenceprocedurestoestimateregressionparameters 14

PAGE 15

andassociatedmeasuresofuncertaintyingeneralizedlinearmodelswithspatiallymisaligneddata.UsingthespatialmisalignmentproblemasamotivatingexampleofBerksonerror,here,weconsidertheBerksonerrorproblemingeneralizedlinearmodelsanddevelopmethodstoestimateregressionparametersandassociatedmeasuresofuncertainty.InChapter 2 ,wereviewandcompareseveralmethodsusedtoconductinferenceinthepoint-to-pointspatialmisalignmentsettingwheretheresponseisassumedtofollowanormaldistribution.InChapter 3 ,weproposeanewmethodtoconductinferenceinspatialmisalignmentproblemswhentheresponsefollowsanormaldistributionandkrigingisusedtoalignthedisparatedatasets.InChapter 4 ,weproposeanewmethodtoconductinferenceaboutregressionparametersingeneralizedlinearmodelswhereerrorinacovariateisdescribedviaaBerksonerrormodel.InChapterChapter 5 wederivetheasymptoticpropertiesofanestimatordevelopedinChapter 3 .InChapter 6 weobtaintheasymptoticpropertiesofaclassofpseudo-likelihoodestimators.Finally,inChapter 7 wedescribeareasforfutureresearch. 15

PAGE 16

CHAPTER2ACOMPARISONOFERRORSINVARIABLESMETHODSFORUSEINREGRESSIONMODELSWITHSPATIALLYMISALIGNEDDATABriefSummary WhenaresponsevariableYismeasuredononesetofpointsandaspatiallyvaryingpredictorvariableXismeasuredonadifferentsetofpoints,XandYhavedifferentsupportsandthusarespatiallymisaligned.TodrawinferenceabouttheassociationbetweenXandY,XiscommonlypredictedatthepointsforwhichYisobserved,andYisregressedonthepredictedX.IfXispredictedusingkrigingorsomeothersmoothingapproach,useofthepredictedinsteadofthetrue(unobserved)Xvaluesintheregressionresultsinunbiasedestimatesoftheregressionparameters.However,thenaivestandarderrorsoftheseparameterstendtobetoosmall.Inthisarticle,twosimulationstudiesareusedtocomparemethodsforprovidingappropriatestandarderrorsinthisspatialsetting.Threeofthemethodsareextendedtothechange-of-supportcasewhereXisobservedatpoints,butYisobservedforarealunits,andtheseapproachesarealsocomparedviasimulation. 2.1BackgroundInmanyenvironmentalepidemiologystudies,thegeographiclocationsoftheexposuremeasurementsdifferfromthelocationsofthehealthdata.Forexample,instudyingthehealtheffectsofairpollution,airqualitymeasurementsarerecordedatmonitorsandhealthinformationisobtainedfromindividualsandassociatedgeographicallywiththeirresidentiallocations.Asanotherexample,healthdataareoftengeographicallyaveragedratesthatmaythenneedtobelinkedwithdatafrom Lopiano,KennethK.,CarolA.Gotway,andLindaJ.Young.2010.Acomparisonoferrorsinvariablesmethodsforuseinregressionmodelswithspatiallymisaligneddata.StatisticalMethodsinMedicalResearch.OnlineFirstasDoi:10.1177/0962280210370266. 16

PAGE 17

environmentalmonitorsatpointlocationsorwithremotely-senseddataatavarietyofresolutions.Suchspatialmisalignmentproblemsareincreasinglyprevalentason-goingprogramsandstudieshaveledtoawealthofinformationthatmaybeusefulforsubsequentenvironmentalepidemiologystudies.Ontheenvironmentalside,theU.S.EnvironmentalProtectionAgency(EPA)collectsairqualitydatathroughanetworkofmonitorsthatcomprisetheAirQualitySystem.TheEPAhasalsocreatedtheSTORETDataWarehouse,arepositoryofwaterqualitymonitoringdatacollectedbynumerouswatermanagementgroups,includingfederalagencies,states,tribes,volunteergroups,anduniversities.TheU.S.DepartmentofAgriculture(USDA)hascollecteddetailedinformationonsoilcharacteristicsthroughouttheUnitedStates.Asforhealthdata,anagencyischargedwithcollectinginformationonhospitaladmissionswithineachstateoftheUnitedStates;forFlorida,theAgencyforHealthCareAdministration(AHCA)hasthisresponsibility.TheU.S.CentersforDiseaseControl(CDC)collectsinformationonavarietyofdiseasesandmaintainsseveraldiseasesurveillancesystems.TheCDC'sBehavioralRiskFactorSurveillanceSystem(BRFSS)providesinformationonavarietyofhealthoutcomemeasuresincludingasthma,obesityandoralhealth,aswellasbehavioralriskfactorsthatmaybeassociatedwiththeseoutcomes.TheNationalCancerInstitute'sSurveillanceEpidemiologyandEndResults(SEER)registryprovidescancerincidenceandsurvivaldatafortheU.S.TheU.S.CensusBureauhason-goingcollectioneffortsonavarietyofsocio-demographicvariablesthatmayalsobeusefulinenvironmentalhealthstudies.Thesearebutafewofthedatasetsthatarepubliclyavailable.Giventheplethoraofavailabledata,programsandstudiesincreasinglyeitheraugmentdatacollectionwithexistingdataorworkentirelywithexistingdatafrommultiplesourcesforanalysisandinference.Becausethepurposesoftheoriginaldatacollectionareoftenquitedifferentfromthoseofthecurrentstudy,datafromdifferentsourcesusuallyhavebeencollectedondifferentgeographicalorspatialunits,andeach 17

PAGE 18

mayvaryfromtheonesofinterest,andthegoalistocombinetheminameaningfulwayforstatisticalinferenceandmodeling.Severalapproachestolinkingsuchdisparatedatahavebeensuggested.Inacase-controlstudyoftheriskofhavingaverylowbirthweightbaby, Rogersetal. ( 2000 )usedanatmospherictransportmodeltoestimateground-levelparticulatematterconcentrationsatresidentiallocationsofmothersidentiedthroughdiseasesurveillanceandvitalstatisticsdatabases.Otherstudieshaveusedsomeformofspatialinterpolationsuchasnearestneighbor(citeMiller2007asanexampleofoneofmanysuchstudies)orkriging( Leemetal. ( 2006 ); Waller&Gotway ( 2004 ))toassignexposuremeasurementstoindividualsortolinkexposuredatawithhealthdata.However,theuseofinterpolatedorpredictedexposuredatainregressionmodelsdesignedtoquantifytheeffectofexposureonhealthpresentsstatisticalchallengesastheuncertaintyassociatedwiththepredictiontechniquesisnotaccountedforwithstandardregressionmethods.Ifapotentialexplanatoryvariablemustbepredictedonadifferentspatialunit,theresultingpredictedvaluesaregenerallysmootherthanthetrueonesand,insomecases,thiscanleadtobiasintheestimatedeffects. Gryparisetal. ( 2009 ); Madsenetal. ( 2008 )Inaddition,smoothingcanfurtherimpactproperuncertaintyassessment.Toovercometheseproblems,severalapproacheswereproposedearlierinthisdecade( Gelfandetal. ( 2001 ); Zhuetal. ( 2003 )),andmorerecentlymeasurementerrormodelshavebeenconsidered( Gryparisetal. ( 2009 ); Madsenetal. ( 2008 ); Szpiroetal. ( 2009 )).Inthestudyconsideredhere,wecompareseveralapproachestoregressionanalysiswithspatiallymisaligneddata.Theseincludethecurrentpracticeofusingpredictedvaluesinregressionanalysisasiftheywereknown,threeadjustedkrigeandregressapproaches( Gryparisetal. ( 2009 ); Madsenetal. ( 2008 ))andadaptationsbasedonbootstrapping( Szpiroetal. ( 2009 )).WeusedataanddesignsimulationstudiesarounddatacollectedinsupportoftheCDC'sEnvironmentalPublicHealth 18

PAGE 19

Tracking(EPHT)program,whichwascreatedtotrackexposureandhealtheffectsthatmayberelatedtoenvironmentalhazards( Young&Gotway ( 2007 ); Youngetal. ( 2009a b )).Thegoalofourcomparisonandevaluationistomakepracticalrecommendationsconcerningthemethodology,concepts,andkeyideassurroundingstatisticalmethodsforspatialmisalignmentproblems.Weinvestigatetheaccuracyofthemethods(bothintermsofunbiasednessandprecisionofestimates),theirsimplicity(intermsofimplementationandcommunication),andtheircomputationalfeasibilityforreal-time,routineusebyepidemiologistsinstatehealthdepartments.Thisistherststudytocomprehensivelyevaluateallofthesemethodsusingthesamesimulationstudyandthersttoevaluatetheminthepoint-to-areachangeofsupportcontext.Although Szpiroetal. ( 2009 )consideredpoint-to-pointmisalignment,theydidnotcomparetheirmethodstothesimplerkrigeandregressmethods,nordidtheyextendtheirmethodstopointtoarealmisalignment.Althoughmoresophisticatedmodelsfortheexposuredatahavebeenconsidered(seeforexample Schoutenetal. ( 1996 )and Sametetal. ( 2000 )),welimitourselvestoanordinarykrigingframework,i.e.allthespatialvariabilityisaccountedforintheerrorcomponentoftheexposurevariable.Asthemethodsweareevaluatinghaveonlybeenimplementedinthecaseofsimplelinearregression,thisisthecontextweconsiderhere.However,describingthehealth-exposurerelationshipthroughrandomeffectsmodelsorspatiallyvaryingcoefcientmodelsismoredesirable.Section 2.2 describesthestatisticalregressionmodelweconsiderhere,aswellasthestatisticalmethodologyweevaluateinthiscontext.Section 2.3 describesthecontextofourworkandthesupportingdata.Section 2.4 givesthedetailsofasimulationstudydesignedtocomparethemethodsdescribedinSection 2.2 .Section 2.5 extendstheproblemtoapoint-to-areamisalignmentscenario,presentstheresultsfortheFloridadata,andpresentstheresultsofadditionalsimulationsinthiscontext.Section 2.6 providesadiscussionofourndingsandplansforfutureresearch. 19

PAGE 20

2.2StatisticalModelsAndMethodsLetYbethevectorofhealthobservationsatlocationswithexposuredata,letXo=(Xo,1(s1),Xo,2(s2),...,Xo,n(sn))0,bethevectorofexposuremeasurements(ozone)observedatpoints(airqualitymonitors).Assumingthatnoneofthehealthobservationsareatlocationswithexposuremeasurements,letXu=(Xu,1(s1),Xu,2(s2),...,Xu,m(sm))0bethevectorofunobservedexposurevaluesatlocationswithhealthobservations.Initially,forsimplicity,weassumeasimplelinearregressionmodelfortherelationshipbetweenenvironmentalexposureandthehealthoutcome: Y=01m1+1Xu+(2)where1m1isam1vectorofones,0and1areunknownregressionparameterstobeestimated,andN(0,)istheerrorofthemodel.Heretheprimaryfocusistheestimationof1anditsstandarderror.Spatialstatisticalmethodsassumethattheenvironmentalexposure(ozone)isarealizationofarandomspatialprocess:fX(s)js2D<2g,whereX(s)isarandomvariableataknownlocations,andsvariessmoothlyoverregionD.AssumethatX(s)hasaconstantmeaninthespatialdomain.Assumethefollowingspatialmodelfortheexposurevariable(ozone): 264Xo(s)Xu(s)375N0B@1n+m1,=264ooouuouu3751CA(2)whereooisthennvariance-covariancematrixofthedata,uuisthemmvariance-covariancematrixamongtheunobservables,andouisthenmvariance-covariancematrixbetweenthedataandtheunobservables.Giventhismodel,thebestlinearunbiasedpredictor(BLUP)oftheunobservablesXuis ^Xu(su)=BLUE()1m1+uo)]TJ /F5 7.97 Tf 6.59 0 Td[(1oo[Xo)]TJ /F1 11.955 Tf 11.96 0 Td[(BLUE()1m1](2) 20

PAGE 21

whereBLUE()isthebestlinearunbiasedestimatoroffromXo, BLUE()=11n)]TJ /F5 7.97 Tf 6.58 0 Td[(1oo 11n)]TJ /F5 7.97 Tf 6.59 0 Td[(1oo1n1Xo.(2)Ingeostatistics,^Xuiscalledthekrigingpredictor( Waller&Gotway ( 2004 )).BecausetheenvironmentalexposureXuin 2 ismorevariablethanitssmoothedpredictor^Xu,wehaveXu=^Xu+v,wherev=Xu)]TJ /F6 11.955 Tf 13.87 2.66 Td[(^Xuistheerrorassociatedwithpredictingexposure.WeassumethatvN(0,v)andisindependentof.Theuseof^XuforXuinmodel 2 resultsinBerksonerror( Gryparisetal. ( 2009 ))asthemodelin 2 canbewrittenas: Y=01m1+1Xu+=01m1+1(^Xu+v)+=01m1+1^Xu+1v+=01m1+1^Xu+(2)where=1v+.Thus,usingthekrigingpredictorof^XuinplaceofXuinducesautocorrelationintotheerrortermofmodel 2 .Intermsofestimationandinference,let^X=(1m1^Xu)and=(0,1).If^XuisanunbiasedpredictorofXu,then ^OLS=(^X0^X))]TJ /F5 7.97 Tf 6.59 0 Td[(1^X0Y,(2)theordinaryleastsquares(OLS)estimatorof,isunbiased.However,asdescribedin Madsenetal. ( 2008 )and Gryparisetal. ( 2009 ),substitutingthepredictedmeanenvironmentalexposuresforthetruevaluesintheregressionexpressionleadstoacorrelatederrorstructure =var()=21v+.(2) 21

PAGE 22

Thus,theOLSvarianceof^1isincorrectbecauseitdoesnotutilizethepropercovariancematrixinEquation 2 thatresultsfromtheadditionaluncertaintyinducedbysmoothing.Howbesttoaccountforthisadditionaluncertainty?Toinvestigatethis,wecomparedthefollowingapproaches: 2.2.1TraditionalKrigeAndRegress(KR)OftenthefactthatXuhasbeenpredictedisignored.Assuming=2I,thepredictedvalues,^Xu,oftheunknownexposurevaluesXuareusedinmodel 2 asiftheywerethetruevalues,resultingintheOLSestimateof1anditsvariancecomputedas^2(^X0^X))]TJ /F5 7.97 Tf 6.59 0 Td[(1[2,2].NoadjustmentfortheadditionaluncertaintyduetospatialpredictionofXuismade. 2.2.2TraditionalKrigeAndRegressWithGeneralCovarianceStructure(KRGC)Usingthepredictedvalues^XuinsteadofthetrueenvironmentalexposurevaluesXuleadstoamorecomplexcovariancestructuregiveninEquation 2 .Perhapsgeostatisticsandvariography,asdescribedin Waller&Gotway ( 2004 ),canbeusedtomodelthespatialstructureintheerrortermandthengeneralizedleastsquarescanbeusedtoestimate1.Inthisapproach ^gls=(^X0^)]TJ /F5 7.97 Tf 6.59 0 Td[(1^X))]TJ /F5 7.97 Tf 6.58 0 Td[(1^X0^)]TJ /F5 7.97 Tf 6.59 0 Td[(1Y,(2)and \var(^gls)=(^X0^)]TJ /F5 7.97 Tf 6.58 0 Td[(1^X))]TJ /F5 7.97 Tf 6.59 0 Td[(1.(2)Thus,thisapproachusesgeneralizedleastsquareswithageneralvariance-covariancematrixinferredthroughgeostatisticaltechniquesandattemptstoaccountforautocorrelationinthedata,boththatarisingnaturallyandthatinducedbysmoothing.Frequently,ageneralcovariancestructureisusedintheregressionprocesstoaccountforanyspatialcorrelationthatremains.Notethat,althoughthemotivationofageneralcovariancestructureisdifferent,themethodisthesamewhethermotivatedfromkriging 22

PAGE 23

theexplanatoryvariableorfrommodelingadditionalspatialcorrelation,whichisnotaccountedforbytheexplanatoryvariables. 2.2.3AdjustedKrigeAndRegressApproachesInadditiontothekrigeandregressapproachesoutlinedabove,weconsideredthetwoadjustedkrigeandregressapproachesin Madsenetal. ( 2008 ) 2.2.3.1Naiveadjustedkrigeandregress(NAKR)InsteadofusingageneralcovariancestructureasinKRGC, Madsenetal. ( 2008 )suggestexploitingtheformofthecovaraincestructurein 2 .Thequasiresiduals,R=Y)]TJ /F6 11.955 Tf 11.96 0 Td[(E(YjXo),havethefollowingcovariancematrix, R=21(u)]TJ /F10 11.955 Tf 11.95 0 Td[(uoo0uo)+.(2)Rcanbeestimatedbyobtaining^R=Y)]TJ /F6 11.955 Tf 13.52 2.59 Td[(^X^andndingtheREMLestimatesofthecovarianceparametersof^R.Thus,isestimatedby ^=^R)]TJ /F6 11.955 Tf 13.43 2.66 Td[(^21,OLS(^u)]TJ /F6 11.955 Tf 13.56 2.66 Td[(^uo^o^0uo).(2)TheKRestimateisnowdenedas ^KR=(^X0^)]TJ /F5 7.97 Tf 6.59 -.01 Td[(1^X))]TJ /F5 7.97 Tf 6.59 -.01 Td[(1^X0^)]TJ /F5 7.97 Tf 6.59 0 Td[(1Y.(2)ThevarianceofthenaiveadjustedKRestimateissimply \var(^KR)naive=(^X0^)]TJ /F5 7.97 Tf 6.59 0 Td[(1^X))]TJ /F5 7.97 Tf 6.59 0 Td[(1.(2) 2.2.3.2Variance-correctedadjustedKR(VCAKR)Thepointestimateisthesameforbothmethodsinthissection,however, Madsenetal. ( 2008 )developamoreaccuratevarianceestimatorasfollows.LetW=[1Xo],where=1mn)]TJ /F10 11.955 Tf 11.95 0 Td[(uo)]TJ /F5 7.97 Tf 6.59 0 Td[(1oo1nn 11n)]TJ /F5 7.97 Tf 6.59 0 Td[(1oo1n1+uo)]TJ /F5 7.97 Tf 6.59 0 Td[(1oo. 23

PAGE 24

Consider ~KR=(W0W))]TJ /F5 7.97 Tf 6.59 0 Td[(1W0Y.(2)Fromthiswehave,~1,KR=X0o0MY X0o0MXowhereM=11m)]TJ /F5 7.97 Tf 6.58 0 Td[(11m1)]TJ /F5 7.97 Tf 6.58 0 Td[(1)]TJ /F10 11.955 Tf 11.96 0 Td[()]TJ /F5 7.97 Tf 6.59 0 Td[(11mm)]TJ /F5 7.97 Tf 6.58 0 Td[(1.Thus, var(~1,KR)=E(Q1Q)]TJ /F5 7.97 Tf 6.59 0 Td[(22)+E(Q3Q)]TJ /F5 7.97 Tf 6.59 0 Td[(22))]TJ /F6 11.955 Tf 11.96 0 Td[((E(Q3Q)]TJ /F5 7.97 Tf 6.59 0 Td[(12))2(2)whereQ1=X0o0M[+21(uu)]TJ /F10 11.955 Tf 11.95 0 Td[(uo)]TJ /F11 7.97 Tf 0 -7.89 Td[(oo10uo)]M0Xo (2)Q2=X0o0MXo (2)Q3=X0o0Muo)]TJ /F5 7.97 Tf 6.59 0 Td[(1ooX0o)]TJ /F7 11.955 Tf 11.96 0 Td[(X0o0Muo)]TJ /F5 7.97 Tf 6.59 0 Td[(1oo1n1 (2)Estimatesoftheseequationscanbefoundbyusingtheestimatesoftheindividualcomponents.TheexpectationscanthenbeestimatedusingMonteCarlointegrationinthefollowingmanner: 1. Repeatedlygeneraterandomvectors^Xofromanormal(BLUE()1n1,^oo) 2. Foreach^Xo,nd^Q1,^Q2and^Q3 3. Calculateaveragesof^Q1^Q)]TJ /F5 7.97 Tf 6.58 0 Td[(22,^Q3^Q)]TJ /F5 7.97 Tf 6.58 0 Td[(22and^Q3^Q)]TJ /F5 7.97 Tf 6.59 0 Td[(12.Thealternativevarianceestimateisthendenedas cvar(^1,KR)=^E(Q1Q)]TJ /F5 7.97 Tf 6.58 0 Td[(22)+^E(Q3Q)]TJ /F5 7.97 Tf 6.59 0 Td[(22))]TJ /F6 11.955 Tf 11.95 0 Td[((^E(Q3Q)]TJ /F5 7.97 Tf 6.58 0 Td[(12))2.(2)However,afterfurtherinspectionweseethatQ2=Q3(seeAppendix).Thisresulthasnon-trivialimplicationsfortheimplementationofthemethodanditisnotevidenttheapproach,initscurrentform,isaviableone. 24

PAGE 25

2.2.4ParametricBootstrap(PB)WiththeKRGCapproach,ageneralcovariancestructureisusedinanattempttoaccountforthemorecomplexerrorstructurearisingfromusingkrigedestimatesinsteadofthetrueX-valueswhenmodelingtheassociationbetweenhealthandenvironmentalexposureasinEquation 2 .However,noconsiderationisgiventotheerrorthatarisesfromusingestimates,insteadofthetrueparametervaluesoftheX-processinthekrigingprocess(seeEquation 2 ).Theerrorthatarisesfromestimatingthekrigingparameters(e.g.,mean,andcovarianceparameters)isclassicalmeasurementerror.UnlikeBerksonerror,classicalmeasurementerrorcancausebiasintheestimateof1aswellasitsstandarderror(see Carrolletal. ( 2006 )). Szpiroetal. ( 2009 )proposedthreebootstrapapproachesthataccountforboththeerrorarisingfromestimatingthekrigingparametersandtheerrorfromusingthesmootherkriged^XvaluesinsteadofthetrueXvaluesinmodelingtheassociationbetweenhealthandexposure 2 .Ineachofthese,ordinaryleastsquaresisusedtoestimate1.Fortherstmethod,thePB,thestandarderrorisfoundbyapproximatingthesamplingdistributionof^1,olsundertheassumedmodelsforenvironmentalexposureandozone.Todothis,bootstrapsamplesaresimulatedusingtheestimatedparametervaluesforenvironmentalexposureandozone,andtheempiricalstandarddeviationof^1,olsfromthesesamplesistakentobethestandarderrorof^1,ols.FormoreconcretenessandassumingMisthenumberofbootstrapsamples,thestepsforthisprocessareasfollowsforagivensetofobservationsYandXo: 1. Estimatetheexposuremodelparameters^and^,inthemodel 2 2. Use^Xuin 2 toestimatemodelparameters^and^. 3. Repeatthestepsbelowforeachj=1,...,M (a) SimulateanewsetofobservationsYjandXo,jbasedonmodelsin 2 and 2 using^,^,^and^. (b) Estimatenewexposureparameters^jand^jusingXo,j. 25

PAGE 26

(c) Derive^Xu,jusing^j,^jandXo,j (d) Calculate^1,jbasedonEquation 2 4. Calculatetheparametricbootstrapstandarderror ^^1=vuut 1 M)]TJ /F6 11.955 Tf 11.96 0 Td[(1MXj=1 ^1,j)]TJ /F6 11.955 Tf 15.76 8.09 Td[(1 MMXj=1^1,j!2.(2) 2.2.5ParameterBootstrap(PRB)Intheparametricbootstrap,theexposuremodelparametersforeachbootstrapsampleareestimatedusingnonlinearoptimization.Inthehopesofreducingcomputationalburden, Szpiroetal. ( 2009 )suggestedtheParameterbootstrap.Heretheexposuremodelparametersforabootstrapsamplearesampledfromtheestimateddensityfunctioncorrespondingtothesamplingdistributionoftheestimatedexposuremodelparameters.Thus,thestepsareasfollows: 1. Estimatetheexposuremodelparameters^and^,inthemodel 2 (a) Estimateadensityfunction^p(,)correspondingtothesamplingdistributionofand^. 2. Use^Xuin 2 toestimatemodelparameters^and^. 3. Repeatthestepsbelowforeachj=1,...,M (a) SimulateanewsetofobservationsYjandXo,jbasedonmodelsin 2 and 2 using^,^,^and^. (b) Sample^and^fromtheprobabilitydistributiondenedby^p(,). (c) Derive^Xu,jusing^j,^jandXo,j (d) Calculate^1,jbasedonEquation 2 4. Calculatetheparameterbootstrapstandarderrorusing 2 2.2.6PartialParametricBootstrap(PPB)InthePB,theexposuremodelparametersforeachbootstrapsampleareestimatedusingnonlinearoptimization.Inthehopesofreducingcomputationalburden, Szpiro 26

PAGE 27

etal. ( 2009 )suggestedthePPB.Heretheerrorinestimatingtheexposuremodelparametersisassumedtobenegligible.Thestepsareasfollows: 1. Estimatetheexposuremodelparameters^and^,inthemodel 2 2. Use^Xuin 2 toestimatemodelparameters^and^. 3. Repeatthestepsbelowforeachj=1,...,M (a) SimulateanewsetofobservationsYjandXo,jbasedonmodelsin 2 and 2 using^,^,^and^. (b) Derive^Xu,jusing^,^andXo,j (c) Calculate^1,jbasedonEquation 2 4. Calculatetheparameterbootstrapstandarderrorusing 2 2.3TheStudyAndSupportingDataTheEPHTprogramisbasedonlittle,ifany,primarydatacollection.Instead,existingenvironmentalexposureandhazardsdata,healthdata,andsocio-demographicdataaretobeusedtoevaluatethespatialandtemporalassociationbetweenhealtheffectsandenvironmentalexposuresandhazards.Linkingdatafromdisparatestudiesonacommonscaleforanalysisandevaluatingtheassociationbetweenhealthandenvironmentalexposuresandhazardsaretwochallengesthattheprogrammustaddress.Forconcreteness,andtoillustrateafewofthemanystatisticalchallengesencounteredinthisbasicstudy,wewillconsidertheassociationbetweentwocoremeasuresoftheEPHTprogram:ozoneandmyocardialinfarction(MI).WorkingwiththedatafromAugust,2005,wewillfocusonmodelingtheassociationbetweenozone(environmentalexposure)andMI(healthoutcome).Ozonemeasurements,recordedfromanetworkofairmonitorsplacedthroughoutthestate(Figure 2-2 ),wereobtainedfromFlorida'sDepartmentofEnvironmentalProtection(FDEP).DuringthestudyperiodofAugust2005,ozonedatawereavailablefrom56monitors.Themaximumofthedailymaximum8-houraverageozonevalues 27

PAGE 28

duringamonthisusedasthemonthlydatavalueforaparticularmonitorandisreferredtoasthemonthlymaximumozone.Florida'sDepartmentofHealth(FDOH)hasadata-sharingagreement(DSA)withFlorida'sAgencyforHealthCareAdministration(AHCA),providinguswithaccesstotwodatasources:condentialhospitalizationrecordsandemergencyroomrecords.ThehospitalizationrecordscontainalladmissionstoFlorida'spublicandprivatehospitalswhereeithertheprimaryorsecondarycauseofadmissionwasMI(ICD-10codes410.0.0)( NationalCenterforHealthStatistics ( 1998 ))Ifapersonpresentstotheemergencyroomandissubsequentlyadmittedtothehospital,thatperson'srecordisincludedinthehospitalizationrecords,butnottheemergencyroomrecords.Non-Floridaresidentswereexcludedfromtheanalysis.ConsistentwithourDSA,AHCAprovidedboththezipcodeandcountyofresidenceforeachpatient'srecord.Selectedpatientdemographicinformationisalsorecorded,includingsex,age,andrace/ethnicity.Foranalysis,thedatamustbelinkedonthesamegeographicalscale,which,forcondentialityreasons,hasbeenthecountylevelinpreviousstudies( Youngetal. ( 2009a b )).TheMIdatawereindirectlystandardizedbyage(aged=45,45,55,and>65years),race/ethnicity(black,white,orother),andsex(femaleormale).Theinformationregardingage,race/ethnicity,andsexfortheMIcaseswereobtainedfromthehospitalrecords.Florida'spopulationwasusedasthecomparisonstandardtocalculatethenumberofexpectedMIcases.ThisgaveanMIstandardizedeventratio(MISER),denedastheratioofthenumberofobservedMIcasestothatexpectedamongtheFloridapopulation,foreachcounty( Woodward ( 2004 )).TolinktheMISERwiththeozonedata,krigingwasusedtopredictozoneontoaregulargridandthenthepredictionswereaveragedwithineachcountytoobtainanaveragemaximumozonevalueforeachcounty.Initially,oursimulationsconsiderthespatialmisalignmentproblembasedonthisgrid(showninFigure 2-1 ).Thepoint-to-areamisalignmentproblemisconsideredinSection 2.5 .Thespatialandstatisticalpropertiesofboththehealth 28

PAGE 29

outcomedataandtheozonedataforthesimulationswerederivedfromthedatainthisexample. 2.4SimulationStudyForMisalignedPointsTwosetsofsimulationswereconductedtocompareandcontrasttheproposedmethodsdiscussedintheprecedingsection.TherstismotivatedbyFlorida'sEPHTprogram.Thesecondismoresimilartothatconsideredby Szpiroetal. ( 2009 ).Inbothcases,thepropertiesoftheestimatorfor1inEquation 2 anditsstandarderrorareofprimaryinterest. 2.4.1Scenario1Forthisrstscenario,itisassumedthatexposureismeasuredatthe56airqualitymonitorsscatteredacrossthestate(Figure 2-2 ).Itisfurtherassumedthathealthvariablesarerecordedatgridpointsacrossthestate(Figure 2-1 ).Foreachof1000simulations,theexposurevariable(ozone)wasgeneratedaccordingtoamultivariatenormaldistributionwithmean42andvariance-covariancematrixdenedbyanexponentialcovariancestructurewithascaleof51andarangeof1,andthetrueXvaluesarerecordedforallmonitorsandgridpoints.ThevaluesofXatthemonitorsareXo,andthoseatthegridpointsareXu.GiventhegeneratedexposurevariableXu,theresponsevariableofinterest,Y,wasgeneratedfollowingthemodelin 2 with0=)]TJ /F6 11.955 Tf 9.3 0 Td[(0.8,1=0.2,and=2.32I.OnceYhasbeengenerated,Xuwasdiscarded,leavingthemisaligneddatasetsXoandY.Eachofthemethodsdescribedintheprevioussectionwereappliedtothe1000datasetstoobtainestimatesof1anditsstandarderror.Theaverage^1andvariance^2^1of^1fromthesimulateddatasets,andtheaverages2^1estimatedvarianceof^1areshowninTable 2-1 .If^1isanunbiasedestimatorof1then^1shouldbewithinerrorofthetruevalueof1=0.2.Ifthevariances2^1of^1isestimatedunbiasedly,s2^1,theaverageoftheestimatedvariancesof^1and^2^1shouldbeapproximatelyequal.ForPB,themethodwasappliedtoonly100ofthe 29

PAGE 30

simulationswithanalM=1000.SincethebootstrapsamplesinthePBcanpotentiallynotconverge,weusedanoriginalMof1250butonlyusedtherst1000thatconverged.ForPPBandPRB,eachmethodwasappliedtoall1000simulationswithM=1000.TheKR,KRGCandNAKRapproacheswereappliedto1000simulations.Allmethodsshowastatisticallysignicantbiasintheestimateofthetrueparameter1=0.2,althoughthebiasineachmethod,withtheexceptionoftheNAKRapproach,maynotbeofpracticalconcernformanyapplications.TheNAKRmethod,anaturalcorrectionbasedontheoreticalderivationsin Madsenetal. ( 2008 ),performedpoorlyinthiscase.Moreover,observingnonpositivedenitecovariancematricesispossiblewiththisapproach.While Madsenetal. ( 2008 )didnotobserveanyintheirsimulations,135ofoursimulationsyieldednonpositive-denitecovariancematrices.Evenafterremovingthese,theobservedcoveragefornominal95%condenceintervalswas51%.Recallthepointestimatesforthenaiveandadvancedmethodsarethesame.Thus,evenaftercorrectingthevarianceasinVCAKR,theperformanceofthesemethodswouldstillbepoor.ThemethodwasoriginallydevelopedforasituationwheretheresponsevariableYconditionalonXisspatiallycorrelated.Sincetheresponsevariableinoursimulationswasnotconditionallyspatiallycorrelated,themethodsintroducedby Madsenetal. ( 2008 )maynotbeappropriate.TheKRandKRGCapproachesunderestimatethevarianceof^1.Asaresult,nominal95%condenceintervalshaveobservedcoverageof90.7%,belowthenominallevel.ThePBandPPBproducedestimatesofthevarianceof^1thatwerebiasedupwards.Consequently,nominal95%condenceintervalshaveanobservedcoverageprobabilityof98.0%and97.1%,respectively.ThePRBmethodconsiderablyoverestimatedtheuncertaintywhencomparedtotheotherbootstrapmethods.Covarianceparameterestimatesinstep3(b)weregeneratedfromatruncatednormaldistributioncenteredattheREMLestimateswithavariance-covariancematrixequaltotheobservedvariance-covariancematrixoftheparameterestimates.Themean 30

PAGE 31

parameterwasgeneratedfromanormaldistributioncenteredatthepointestimatewithvarianceequaltotheestimatedvarianceoftheparameterestimate.ThePRB,althoughessentiallyunbiased,overestimatesthevariance.Itislikelythatthedisconnectcreatedwhensimulatingdatasetsbasedonacertainsetofparametervaluesandthenrandomlygeneratingparametervaluestotthedatasetsisdrivingtheseinatedestimates.Basedonthislimitedstudy,thePPBappearstoperformthebest. 2.4.2Scenario2Forthesecondsimulationstudy,a400500gridwasestablished.Fromthegrid,200pointswererandomlyselectedtobepointsforwhichexposuredatawererecorded,and300pointswererandomlyselectedaspointsforwhichhealtheffectswereobserved.Thedistributionusedtogeneratetheexposurevariableisamultivariatenormaldistributionwithmean1.36andvariance-covariancematrixdenedbyanexponentialcovariancestructurewithascaleof3.76,arangeof24.13andanuggetof1.34.Giventhegeneratedexposurevariable,theresponsevariableofinterestYwasgeneratedfollowing 2 with0=5.061=)]TJ /F6 11.955 Tf 9.3 0 Td[(0.322and=0.76I.Thesimulationstudyconsistedofgenerating2500exposuresurfaces.Ofthese2500surfaces,onlyasubsetwasselected.Thissubsetwaschosenbasedonthefollowingrulestoensuretheaccuracyoftheinitialtsoftheexposuresurface: IfanydiagonalelementoftheHessianmatrixwasgreaterthan9thendelete; Iftheestimateofthenuggeteffectislessthan0.05thendelete; Iftheestimateoftherangeparameterislessthan5orgreaterthan200thendelete; Iftheestimateofthescaleparameterisgreaterthan7thendelete.Theparametervalueswerechosentocoincidewiththoseusedin Szpiroetal. ( 2009 )Similarly,therulesfordeletionwerebasedontherecommendationsof Szpiroetal. ( 2009 ) 31

PAGE 32

ForthePB,themethodwasappliedtoonly100ofthesimulationsandanalM=1000.Bootstrapsamplesthatdidnotconvergewerediscarded.Also,severalofthebootstrapsamplesconvergedpoorlyandyieldedlargeestimatesforthevarianceof^1,j.Whenthevarianceof^1,jwasgreaterthan1,thatbootstrapsamplewasremoved.ForPPBandPRB,themethodwasappliedto2000simulationswithM=1000.TheKRGCdidnotalwaysconverge;thus,theresultsforthismethodarebasedonthe1714simulationsthatdidconverge.TheKRandNAKRprocedureswereappliedto2000simulations.TheresultsforeachofthemethodsappliedinthissimulationscenarioareshowninTable 2-2 .Theresultsinthissimulationstudyaresimilartothosepresentedintheprevioussection.Allmethodsexhibitedaslightbiasintheestimateofthetrueparameter1=)]TJ /F6 11.955 Tf 9.29 0 Td[(0.322.TheKRandKRGCapproachesagainunderestimatethevarianceof^1.Asaresult,nominal95%condenceintervalshaveobservedcoverageof89.1%and91.6%respectively,wellbelowthenominallevel.ThePBandPPBhadestimatedvariancesof^1thatslightlyoverestimatethetruevariance.Consequently,nominal95%condenceintervalshaveanobservedcoverageof98.0%and95.8%respectively.ThePRBmethodconsiderablyoverestimatedtheuncertaintywhencomparedtotheotherbootstrapmethods.Covarianceparameterestimatesinstep3(b)weregeneratedfromatruncatednormaldistributioncenteredattheREMLestimateswithavariance-covariancematrixequaltotheobservedvariance-covariancematrixoftheparameterestimates.Themeanparameterwasgeneratedfromanormaldistributioncenteredatthepointestimatewithvarianceequaltotheestimatedvarianceoftheparameterestimate.ThePRB,althoughessentiallyunbiased,overestimatesthevariance.Thereasonsforthispoorperformancenotedintheprevioussectionarelikelythereasonforpoorperformanceinthisscenarioaswell.TheNAKRmethodperformedpoorlyinthiscaseaswell.Havingremoved246simulationsthatyieldednonpositive-denitecovariancematrices,nominal95% 32

PAGE 33

condenceintervalshaveanobservedcoverageof51.4%.Thereasonsforthepoorperformanceinthepreviousscenarioarelikelyreasonsforpoorperformanceinthisscenarioaswell.Thus,themethodsintroducedby Madsenetal. ( 2008 )arenotviablefortheapplicationathand.Again,itappearsthatthePPBmethodperformsthebesthere. 2.5MovingFromPointstoAreasFortheFloridastudy,healthdataareavailableattheareal(county)level,notthepoint(residential)level.Thus,thesimulationstudyofmisalignedpointdataprovidesonlysomeinsightintoaproblemthatoccurswithmanyEPHTstudies.Inwhatfollows,weconsideruncertaintyassessmentinaregressionmodelwhereexposureismeasuredatpointlocationsbutthehealthdataareratesassociatedwithgeographicareas(counties).WebeginrstbyapplyingeachofthemethodsdescribedinSection 2.2 totheFloridadata,andthenfurtherassesstheaccuracyofthesemethodsviatwosimulationstudies.Therstofthesebuildsonscenario1ofthemisalignedpointsimulation.Thesecondbuildsonscenario2inSection 2.2 2.5.1ResultsForTheFloridaDataBecauseenvironmentalexposureismeasuredonthepointscale,butinferenceistobeconductedforarealunits,aspatialpredictiontechniqueknownasblockkrigingisusedtopredictenvironmentalexposureforthearealunits.Withthistechnique,krigingwasusedtopredictozoneontoaregulargridasdescribedpreviously,butthenthepredictionswereaveragedwithineachcountytoobtainanaveragemaximumozonevalueforeachcounty.Giventhepredictedarealaveragesoftheozonedata,thehealthoutcomedataareregressedonthepredictedenvironmentalexposureasinEquation 2 .Theregressionsarethusbasedoncounty-levelinformation.FiveofthemethodsinSection 2.2 thatcanbeeasilyextendedtopointtoarealsettingswereusedtoanalyzetheFloridadatadiscussedinSection 2.3 (seeTable 33

PAGE 34

2-3 ).A12-kmgrid,andnotthecoarsergridfromthesimulation,wasusedforourinitialinvestigationasitgreatlyreducedcomputationaltime.Allthemethodsyieldedsimilarpointestimatesfortheparameterofinterest,1.Thestandarderrors,however,aredifferentforeachofthemethods.TheKR,KRGC,PRBandPPBallhavesimilarestimatesofthestandarderror,althoughthestandarderrorobtainedusingthePBwasalmosttwiceaslarge.Toinvestigatethisfurther,additionalsimulationswereconducted.Asbefore,weconsidertwoscenariosdescribedbelow.However,onlythreeofthemethods(KR,KRGC,andPPB)wereconsideredsincethePRBandPPBapproacheswerecomputationallyintensivetoimplementinthecomplexsimulationstudyconsideredhere. 2.5.2FloridaGridSimulationTheexposuredataweregeneratedforthemonitorsandgridpointsasforscenario1inSection 2.4 .Theaverageoftheexposurevaluesatgridpointswithinacountywasusedtoobtainx(Bi),theaverageexposureforthecounty.Then,conditionalonthiscounty-averagedexposure,theresponsevectorY(B)wasgeneratedaccordingtoamultivariatenormaldistributionwithmeanwith)]TJ /F6 11.955 Tf 9.3 0 Td[(0.8+0.2x(B)andvariance-covariancematrix=2.32I.Theexposurevaluesatthegridpointswerediscarded,andblockkrigingwasusedtopredictenvironmentalexposureattheareallevel,^xB.Then,theresponseY(B)wasregressedon^xB.Theregressionwasconductedintwodifferentways,onewithanindependenterrorstructure(KR)andonewithageneralcovariancestructurebasedonanexponentialcovariancestructure(KRGC).WealsoassessedtheperformanceofPPBextendedtoapointtoarealsetting.Theresultsbasedon1000simulationsaregiveninTable 2-4 .ThePPBwasimplementedusinganM=1000.Itisinterestingtonotethebiasobservedintheparameterestimateof1isnotofpracticalconcern.Thevarianceestimatesseemfairlyaccurate,althoughthecoverageprobabilitiesarenotentirelyaccurateforanyofthemethods.Nominal95%condenceintervalsfortheKRandKRGCmethodshaveanobservedcoverageof93.6%and 34

PAGE 35

93.7%respectively.ThePPBoverestimatesthevariability.Asaresult,nominal95%condenceintervalshaveobservedcoverageof97.3%. 2.5.3Scenario2SimulationFirst,the400500gridwaspartitionedinto32arealunits.Usingtheacceptedexposuredatasets,averageswithineacharealunit,x(Bi),i=1,...,32wereobtainedandtheresponsevectorY(B)wasgeneratedfollowingamultivariatenormaldistributionwithmean5.06)]TJ /F6 11.955 Tf 13.49 0 Td[(0.322x(B)andvariance-covariancematrix=0.76I.Theexposurevaluesatthegridpointswerediscarded,andblockkrigingwasusedtopredictenvironmentalexposureattheareallevel,^xB.Then,theresponseY(B)wasregressedon^xB.Theregressionwasdonethreedifferentways,onewithanindependenterrorstructure(KR),onewithageneralcovariancestructurebasedonanexponentialcovariancestructure(KRGC)andonebasedonextendingthepartialparametricbootstrap(PPB).Theresultsbasedon2000simulationsaregiveninTable 2-5 .ThePPBwasimplementedusinganM=1000.Theparameterestimatesandvarianceestimatesappeartoperformquitewellinthisscenario.Thebiasintheparameterestimateandinitsvarianceestimatearenotofpracticalimportanceineithercase.Consequently,thenominal95%condenceintervalshaveanobservedcoverageprobabilityofabout95%inbothkrigeandregresscasesand95.9%inforthePPB.Inthisscenario,itisclearthattheKRandKRGCmethodsperformadequatelyandtheextraeffortinthePPBisnotnecessary.Itisinterestingtonotethatthestandarderrorsoftheregressioncoefcientislargerwhenthemisalignmentoccursfrompointstoarea. 2.6ConcludingRemarksIncreasingly,publichealthisassociatedwithenvironmentalexposures.Becausethepublichealthandenvironmentaldataareseldomcollectedatthesamepointsintimeandspace,issuesofchange-of-supportarecommon( Gelfandetal. ( 2001 ); Mugglin&Carlin ( 1998 ); Mugglinetal. ( 1999 )).Priorto Madsenetal. ( 2008 )and Gryparisetal. 35

PAGE 36

( 2009 ),noconsiderationhadbeengiventothepotentialbiasthatcouldbeintroducedusingpredictionsatpointsinsubsequentspatialmodels.Themainobjectivesofthispaperwereto:1)evaluatetheexistingerrorsinvariablesmethodologytoaddresspoint-to-pointspatialmisalignmentproblems;and2)understandthebehaviorofparameterestimatesobtainedfromtheseapproachesinapoint-to-areamisalignmentsetting.Fortherstobjective,oursimulationstudiesindicatethatallapproachesgiveanessentiallyunbiasedestimateoftheregressionparameter,1,withtheexceptionoftheNAKR.Thenaiveadjustedkrigeandregressmethodandthevariancecorrectedadjustedkrigeandregressmethoddevelopedby Madsenetal. ( 2008 )donotappeartobeviablesolutions.Theotherkrige-and-regressapproachesdonotfullyaccountforthevariabilityinpredictedexposuredataresultingincoverageprobabilitiesthataretoolow.Thebootstrappingapproachesseemtoover-adjustforthisvariability,producingcoverageprobabilitiesthataretoolarge.SomeattentioncouldbegiventodevelopingawaytosystematicallyaccountforpoortsinthePRBbootstrapsamplestoimprovetheperformanceofthemethod,butitisnotcleariftherulesforonesetofdatawillbethesameforanothersetofdata.Similarly,awaytosystematicallyremovebootstrapsamplesinthePBisdesirable.Itisnotimmediatelyclearhowtodealwithconvergenceissuesamongthebootstrapsamplesingeneral.Thisproblemmayhinderthemethods'utilityincertainapplications.HowtohandlethisproblemingeneralisperhapsanavenueoffutureresearchshouldthePBbeadoptedforanalysis.Decidingonawaytoremovesamplesmayvaryacrossapplications.Forthisreason,usingthePRBorPBinanautomatedfashionmightnotberealistic.Also,removingsamplesmayintroduceabiasinthevarianceestimate.ThePPBmaybethemostaccurateofallapproachesconsideredhereandappearstowarranttheadditionalcomputationaleffortwhenappliedtomisalignedpointdata.Becausetheadditionalcomputationaleffortofapplyingthemethodspresentedin Szpiroetal. ( 2009 )toonedatasetisnegligible,ourmainconcernwiththesemethodologiesisnotthe 36

PAGE 37

computationaltime,butinsteadtheperformanceofthemethods.Basedonthislimitedstudy,itisevidenttheperformanceofthemethodsappearstodependonthesimulationscenario.Withrespecttooursecondobjective,theresultswereparticularlysurprising.Weexpectedunbiasednessoftheparameterestimate,however,wedidnotexpecttoseeunbiasednessintheestimateofthestandarderror.Apparentlytheaveraging(andperhapstheautocorrelationinducedbyit)increasesthestandarderrorsoftheregressioncoefcient.Thus,itisclearincertainsituationswithpointtoarealmisalignment,thebasicKRmethodperformswell.However,itisnotclearthiswillbetrueingeneral.GriddensityandotherfactorsthatcontributetoobtainingblockkrigingestimatesmayaffecttheperformanceoftheKRmethod.Futureresearchwilllikelyshedmorelightontheperformanceoftheseapproachesinpoint-to-areamisalignmentproblems.Thereareseveralavenuesforfutureresearch.First,allofthesemethodswereevaluatedassumingaconstantmeanfortheexposuresurface.Assumingthatallofthespatialvariabilityisaccountedforbythespatialautocorrelation,i.e.assumingaconstantmeanmodel,presentsaninterestinglimitation.Iftheexposurevariableisconstantacrosstheregionofinterest,thenthereisnowaytoestimatetheeffectofexposureontheresponsevariable.Therefore,thespatialvariabilitymustcreateenoughvariabilityintheexposurevariabletorealisticallyestimatetheeffecttheresponsevariable.Ifatrendispresentoverthestudyarea,thentheexposuresaremorelikelytohavearangethatallowsforbetterestimationoftheassociation.However,aslighttrendsurfacewithlittletonospatialvariationcouldbejustasproblematic.Nonetheless,itisnotclearhowtheresultswoulddifferifamorecomplicatedmeanfunctionwereconsidered,andthebootstrappingmethodsmayneedtobeadaptedtoincludethisscenario.Somethoughtonhowbesttoadaptalltheerrorsinvariablesapproachestochangeofsupportproblemsotherthanthepoint-to-pointmisalignmentproblemisneeded. 37

PAGE 38

Finally,comparisontoBayesianhierarchicalmodelsthatarephilosophicallyappealinginthesemodelingscenariosisneeded.ABayesianapproachtospatialmisalignment,seeforexample Banerjee&Gelfand ( 2002 ),isofparticularinterestsinceitnaturallyaccountsfortheextrauncertainty.TheBayesianapproachappearstoworkrelativelywellinthescenariosin Gryparisetal. ( 2009 ).However,weencounteredissueswhenimplementingaBayesianmodelusingcommonsoftware,namelyJAGSandWinBUGS.AlthoughtheobvioussolutionistodevelopaGibbssampler,most,ifnotall,statehealthdepartmentslackthetechnicalcapacitytodothis.DevelopingaGibbssamplerorwritingWinBUGSscriptsforeveryparticularapplicationisachallengingtask,evenforstatisticians.Thus,aBayesianapproachisnotaviableoptionforroutineapplieduse.Nonetheless,anevaluationofBayesianhierarchicalmodelsisneeded,althoughthisgreatlyincreasesthecomputationalburdenofalreadydemandingsimulationstudies. 2.7SupportWearegratefultoDr.LisaMadsenandDr.AdamSzpiroforansweringquestionsregardingtheirmethodology.WealsoappreciateDr.GregKearneyandChrisDuclosattheFloridaDepartmentofHealthforprovidingtheirsupportofourresearch.ThispublicationwassupportedbytheFloridaDepartmentofHealth,DivisionofEnvironmentalHealthandGrant/CooperativeAgreementNumber5U38EH000177-02fromtheCentersforDiseaseControlandPrevention(CDC).ThendingsandconclusionsinthisreportarethoseoftheauthorsanddonotnecessarilyrepresenttheviewsoftheCentersforDiseaseControlandPrevention. 38

PAGE 39

Table2-1. Scenario1SimulationResults:Average^1andVariance^2^1oftheEstimatesof1andAverageoftheEstimatedVariances2^1of^1fortheSimulatedDatasetsinScenario1 EstimationMethod^1^2^1s2^1CoveragePercentage KR0.2080.003730.0026590.7KRGC0.2080.003730.0026690.7NAKR0.1702.1170.0037151.0PB0.2080.005050.0055298.0PRB0.2080.003730.97899.7PPB0.2080.003730.0042397.1 Table2-2. Scenario2SimulationResults:Average^1andVariance^2^1oftheEstimatesof1andAverageoftheEstimatedVariances2^1of^1fortheSimulatedDatasetsinScenario2 EstimationMethod^1^2^1s2^1CoveragePercentage KR-0.3260.005750.0037389.1KRGC-0.3290.006840.0045791.6NAKR-0.33012.5410.0029351.4PB-0.3200.004670.0083798.0PRB-0.3260.0057561.58999.0PPB-0.3260.005750.0062995.8 Table2-3. Estimatesof^1andtheStandardErrorof^1forFlorida EstimationMethod^1s^1 KR0.0240.018KRGC0.0280.020PB0.0240.048PRB0.0240.025PPB0.0240.017 39

PAGE 40

Table2-4. Scenario2SimulationResults:Average^1andVariance^2^1oftheEstimatesof1andAverageoftheEstimatedVariances2^1of^1fortheSimulatedDatasetsintheScenario2point-to-arealSimulation EstimationMethod^1^2^1s2^1CoveragePercentage KR-0.3250.03230.033595.0KRGC-0.3250.03230.033595.0PPB-0.3250.03230.036895.9 Table2-5. Average^1andVariance^2^1oftheEstimatesof1andAverageoftheEstimatedVariances2^1of^1fortheSimulatedDatasetsintheFloridaGridSimulation EstimationMethod^1^2^1s2^1CoveragePercentageKR0.2040.00530.005093.6KRGC0.2050.00530.005193.7PPB0.2040.00530.006797.3 Figure2-1. CoarsegridlaidoverthestateofFlorida 40

PAGE 41

Figure2-2. AirmonitorlocationsinthestateofFlorida 41

PAGE 42

CHAPTER3SPATIALCOVARIANCEESTIMATIONINSPATIALLY-MISALIGNEDREGRESSIONMODELSWITHBERKSONERRORBriefSummary Inenvironmentalstudies,relationshipsamongvariablesthataremisalignedinspaceareroutinelyassessed.Becausethedataaremisaligned,krigingisoftenusedtopredictthecovariateatthelocationswheretheresponseisobserved.UsingkrigingpredictionstoestimateregressionparametersinlinearregressionmodelsintroducesBerksonerror.Asaresult,theBerksonerrorinducesacovariancestructurethatischallengingtoestimate.Inaddition,iftheparametersassociatedwiththecovariateareestimated,thenanaspectofclassicalmeasurementerrorisintroduced.WecharacterizethemeasurementerroraspartofabroaderclassofBerksonerrormodelsanddevelopanestimatedgeneralizedleastsquaresestimatorusingestimatedcovarianceparameters.Inworkingwiththeinducedmodel,wefullyaccountfortheerrorstructureandestimatethecovarianceparametersusinglikelihoodbasedmethodsandmethodofmoments.WeassesstheperformanceoftheestimatorsusingsimulationandillustratethemethodologyusingpubliclyavailabledatafromtheU.S.EnvironmentalProtectionAgency.Finally,weextendtheresultstoanotherchange-of-supportproblemwheretheresponseisobservedatthearealunitlevelandthecovariateisobservedatthepointlevelusinganexamplefromtheCentersforDiseaseControl'sEnvironmentalPublicHealthTrackingProgram. 3.1BackgroundEnvironmentalepidemiologistsassessrelationshipsbetweenenvironmentalexposureandhumanhealthoutcomes( Gryparisetal. ( 2009 )).Ecologistsstudytherelationshipsamongmeasuresofecosystemqualityandspatiallyvaryingcovariates( Zhangetal. ( 2011 )).Moregenerally,therelationshipsamongvariableswithanexplicitspatialcomponentareroutinelyassessed.Regardlessoftheapplication,thevariablesareoftenobservedatdifferentlocationsoraggregatedoverdifferent 42

PAGE 43

geographicalunits.Suchdataaresaidtobespatiallymisaligned.Forexample,apoint-to-pointmisalignmentproblemariseswhenrelatingairqualitymeasurements,observedatmonitors(points),andbirthweightsobservedattheresidentiallocationsofthemothers(differentpoints).Asanotherillustration,apoint-to-arealmisalignmentproblemariseswhenrelatingozonemeasurements,observedatmonitors(points),andratesofmyocardialinfarction,observedatthecountylevel(arealunits).Here,thegeneralpoint-to-pointmisalignmentproblemwheretheresponsevariableisobservedatonesetofpointsinspaceandtheexplanatoryvariableisobservedatanothersetofpointsinspaceisconsideredrst.Then,theresultsareextendedtothepoint-to-arealmisalignmentproblemwheretheresponsevariableisobservedonasetofarealunitsandtheexplanatoryvariableisobservedatasetofpoints.SupposetheresponsevariableYisobservedatnspatiallocationssu,andtheexplanatoryvariableXisobservedatnodifferentspatiallocationsso.Initially,considerasimplelinearregressionmodelfortherelationshipbetweenYandX: Y(su)=01n1+1X(su)+,N(0,()),(3)whereX(su)istheunobservedexplanatoryvariableatlocationssuandisavectorofparametersthatdenethestructureof.Inpractice,predictedvaluesofX(su)areusedtoalignthedisparatedatasetsandestimatetheparametersin 3 .Supposetherelationshipbetweenthepredictedvalues~X(su)andthetruevaluescanbewrittenas~X(su)=X(su)+x,wherexN(0,x).If~X(su)isusedasthecovariatein 3 ,thenaclassicalmeasurementerrormodelresultsandordinaryleastsquares(OLS)estimatesof1arebiased.See Carrolletal. ( 2006 )or Buonaccorsi ( 2009 )forareviewofmethodsthataccountforclassicalmeasurementerror.Althoughimportant,classicalmeasurementerrormethodsarenotthefocushere.Inpractice,smoothingtechniques,suchaskrigingandgeneralizedadditivemodels,arefrequentlyusedtoaligndisparatedatasets( Gryparisetal. ( 2009 ); Pacioreketal. 43

PAGE 44

( 2009 )).Asnotedinearlierwork,seeforexample Pacioreketal. ( 2009 ),smoothedestimates,SX(su),canberelatedtothetruevaluesthroughaBerksonerrorrelationship: X(su)=SX(su)+v,(3)wherevisameanzerorandomvariablethatisorthogonaltoSX(su)(e.g.,vN(0,v).Althoughintheclassicalmeasurementerrormodelthepredictedvaluesaremorevariablethanthetruevalues,intheBerksonerrormodelthesmoothedestimatesarelessvariablethantheactualvalues.Becauseoftherelationshipin 3 Y(su)=01+1X(su)+=01+1(SX(su)+v)+=01+1SX(su)+1v+. (3) Asaresult,ifthesmoothedvaluesdependonasetofparametersandtheparametersareknown,thentheBerksonerrormodeldoesnotcauseOLSestimatesof1tobebiased.However,ifthesmoothedvaluesareestimatedusinganestimateof,thenclassicalmeasurementerrorarisesandcancauseOLSestimatestobebiased.Moreover,theadditionalerrorterm1vresultsinacorrelatederrorstructurethatneedstobeaccountedforinestimatesofuncertainty.Inthispaper,wedevelopacomputationallyefcientestimatedgeneralizedleastsquaresestimatorof1intheinducedmodelthataccountsfortheerrorstructureinducedbythepredictionofmisaligneddata.Morespecically,weconsiderthecasewhenkrigingpredictionsareusedtoalignthedisparatedatasets.Thetwo-stageapproachwedevelopinSection 3.3 presentsacomputationallyefcientalternativetofullyBayesianmethods( Gryparisetal. ( 2009 ); Mugglinetal. ( 2000 )).AlthoughothershaveconsideredGLSestimationinthiscontext,theestimatedcovariancematricesneededforroutineimplementationhavebeenbasedonapproximationsthatmayperformpoorlyinpracticalapplications.Morespecically, Madsenetal. 44

PAGE 45

( 2008 )introducedageneralizedleastsquaresestimatorbasedontheoriginalmodelinEquation 3 .AsshowninSection 3.3 ,theproposedvarianceoftheestimatorcanperformpoorlydependingontheinducederrorstructure.Inaddition, Gryparisetal. ( 2009 )discussageneralizedleastsquaresestimatorintheonlinesupplementaryleusingpenalizedregressionsplinestopredictthecovariate.AsdiscussedinSection 3.3 ,thisapproachusuallyresultsinclassicalmeasurementerroraswellasBerksonerror.Finally, Szpiroetal. ( 2011 )partitionedtheerrorintoaclassical-likecomponentandaBerkson-likecomponentandintroducedthreebootstrapmethodologiestoaccountforbothsourcesoferrorinanordinaryleastsquaresestimator.Thecovarianceparameterswereestimatedusingtheordinaryleastsquaresttoanempiricalbinnedsemivariogramoftheresidualsfrom 3 .Becausetheseresidualsareassociatedwiththeerrorterm1v+,inSection 3.3 weshowthetheoreticalsemivariogramisafunctionofthedistancebetweentwolocationsandthedistancebetweenthetwolocationsandthelocationsso.Therefore,theempiricalbinnedsemivariogramasafunctionofdistancealonedoesnotaccuratelyreecttheinducederrorstructure.Inlightofthepreviousworkdiscussedhere,aneedremainsinspatialmisalignmentproblemsformethodstoestimateand,methodstoincorporateestimatedcovarianceparametersintoavalid,computationallyefcientEGLSestimatorof.WecasttheprobleminthecontextofabroaderclassofBerksonerrorproblemsandaddressbothneedsusingatwo-stageprocedure.Intherststage,theparametersassociatedwithXareestimatedbasedonX(so).Inthesecondstage,theseestimatedparametersareassumedknownandisestimatedusingmethodofmomentswhen=2Iandmaximumlikelihoodappliedtoerrorcontrastswhenhasamorecomplicatederrorstructure.Then,anestimatedgeneralizedleastsquaresestimatorofisobtainedbasedontheinducedmodel 3 .Moreover,weextendtheideastothepoint-to-arealchangeofsupportproblemmorecommonlyencounteredinhealthapplicationswherehealthdataareaggregatedovergeographicunits.Thepropertiesoftheapproach 45

PAGE 46

areassessedusingasimulationstudyandillustratedintwoexamplesusingpubliclyavailabledatafromtheU.S.EnvironmentalProtectionAgency'sEnvironmentalMonitoringandAssessmentProgramanddatafromtheCentersforDiseaseControl'sEnvironmentalPublicHealthTrackingProgram.Thepaperisorganizedasfollows.InSection 3.2 ,weintroducethenotationandmodelingframework.TheproposedmethodologyispresentedinSection 3.3 andappliedtotwodatasets,onefromtheU.S.EPAandtheotherfromtheCDC,inSection 3.4 .InSection 3.5 ,wepresentasimulationstudy.Finally,weconcludewithadiscussioninSection 3.6 3.2ModelingFrameworkTomakethemodelingframeworkmoregeneral,additionalcovariatesobservedatthesamepointsastheresponseareconsideredinmodel 3 andtheproblemisrecastinamultiplelinearregressionframework.LetZ(su)beanncmatrixwhereeachrowcorrespondstoccovariatesobservedateachofthenlocationsinsu.SupposetherelationshipamongY(su),X(su),andZ(su)is Y(su)=01n1+1X(su)+Z(su)[2,...,c+1]0+,(3)whereN(0,())andisaanunknownvectorofmparametersthatdenethecovariancematrix.Inaddition,supposeX(su)andX(so)arejointlydistributedas 264X(su)X(so)375N0B@264C(su)C(so)375,x(x)=264uu(x)uo(x)ou(x)oo(x)3751CA,(3)whereC(su)andC(so)arenkandnokmatricesofobservedcovariatesatthelocationsinsuandsorespectively,isanunknownvectorofkparameters,xisanunknownvectoroflparametersthatdenethecovariancematrixuuamongtheunobservedlocations,thecovariancematrixuoamongtheunobservedandobserved 46

PAGE 47

locations,andthecovariancematrixooamongtheobservedlocations.Forexample,xmaybeparametersthatdeneanexponentialcovariancestructure.Theconditional(kriging)distributionofX(su)jX(so)is X(su)jX(so),,xN(X(su)=C(su)+uo)]TJ /F5 7.97 Tf 6.59 0 Td[(1oo(X(so))]TJ /F10 11.955 Tf 11.95 0 Td[(C(so)),=uu)]TJ /F10 11.955 Tf 11.96 0 Td[(uo)]TJ /F5 7.97 Tf 6.58 0 Td[(1ooou). (3) NoticethekrigingmeanhasaBerksonerrorrelationshipwiththeunobservedcovariate, X(su)jX(so),,x=X(su)+v,(3)wherevN(0,).Asaresult,themodel Y(su)jX(su),Z,,=01n1+1X(su)+Z(su)[2,...,c+1]0+ (3) X(su)jX(so),,x=X(su)+v (3) followsaBerksonerrormodel.Here,followingtheterminologyof Buonaccorsi ( 2009 ),theBerksonerrormodelisdenedbythesetofBerksonparameters. def=f,xg.(3) 3.3ProposedMethodologyFormeasurementerrorproblems,itisusefultoconsidertheinducedmodel Y(su)jX(so),Z,,,=0+1X(su)jX(so),+Z(su)[2,...,c+1]0+=0+1(X(su)+v)+Z(su)[2,...,c+1]0+=0+1X(su)+Z(su)[2,...,c+1]0+1v+. (3) Themeanandvarianceoftheinducedmodelare E(Y(su)jX(so),Z,,,)=01n1+1X(su)+Z(su)[2,...,c+1]0def=X, (3) 47

PAGE 48

whereX=[1n1X(su)Z(su)]and=[012...c+1]0,and Var(Y(su)jX(so),Z,,,)=21Var(v)+Var()=21+def=y. (3) BecauseX(su)jX(so),andarenormallydistributed,Y(su)jX(so),Z,,,isnormallydistributed.Thatis, Y(su)jX(so),Z,,,=X+,(3)whereN(0,y).Thegeneralizedleastsquaresestimator,andhencethebestlinearunbiasedestimator,ofis ^GLS=(X0)]TJ /F5 7.97 Tf 6.59 0 Td[(1yX))]TJ /F5 7.97 Tf 6.58 0 Td[(1X0)]TJ /F5 7.97 Tf 6.59 0 Td[(1yY(su).(3)Thevarianceofthisestimatoris Var(^GLS)=(X0)]TJ /F5 7.97 Tf 6.59 0 Td[(1yX))]TJ /F5 7.97 Tf 6.58 0 Td[(1.(3)Asthenaturalchoiceofanestimatorinthismodel,theGLSestimatorhasbeenidentiedandbrieydiscussedinpreviouswork( Gryparisetal. ( 2009 ); Szpiroetal. ( 2011 )).However,explicitdetailsregardingtheGLSestimatoranditsimplementationinpracticehavenotyetbeenpresented.InadditiontotheGLSestimator,severalweightedleastsquaresestimatorshavebeendiscussed( Gryparisetal. ( 2009 ); Madsenetal. ( 2008 )).Notably,theweightedleastsquaresestimator ^WLS=(X0)]TJ /F5 7.97 Tf 6.59 0 Td[(1X))]TJ /F5 7.97 Tf 6.58 0 Td[(1X0)]TJ /F5 7.97 Tf 6.59 0 Td[(1Y(su)(3) 48

PAGE 49

wasconsideredpreviously( Madsenetal. ( 2008 )).Althoughanassociatednaivevarianceestimatordenedas Var(^WLS)=(X0)]TJ /F5 7.97 Tf 6.59 0 Td[(1X))]TJ /F5 7.97 Tf 6.58 0 Td[(1(3)hasperformedwellinsimulationstudies( Madsenetal. ( 2008 )),because Var(^WLS)=(X0)]TJ /F5 7.97 Tf 6.58 0 Td[(1X))]TJ /F5 7.97 Tf 6.59 0 Td[(1X0)]TJ /F5 7.97 Tf 6.58 0 Td[(1y)]TJ /F5 7.97 Tf 6.59 0 Td[(1X(X0)]TJ /F5 7.97 Tf 6.59 0 Td[(1X))]TJ /F5 7.97 Tf 6.58 0 Td[(1 (3) 6=(X0)]TJ /F5 7.97 Tf 6.58 0 Td[(1X))]TJ /F5 7.97 Tf 6.59 0 Td[(1, (3) thevariancewillnotbecorrectingeneral.If)]TJ /F5 7.97 Tf 6.59 0 Td[(1y=I,thenthevarianceestimatorshouldperformwellinestimatingtheassociateduncertainty.However,thiswillnotbetrueingeneral. 3.3.1GLSEstimationFornow,theBerksonparametersareassumedknownandmethodstoconductinferencefor1in 3 aredevelopedaccordingly.Then,inferenceproceduresforwhentheparametersareunknownaredeveloped.Toestimateasin 3 ,anestimateofymustbefound.Thisischallengingbecause,asseenin 3 ,ydependsnotonlyonandbutalsoon1,theunknownparameterofinterest.Here,weworkexplicitlywiththekrigingdistributionusingthemean(i.e.,thebestlinearunbiasedpredictor(BLUP))andvariancein 3 resultingintheinducedmodelin 3 .Asshownin 3 ,aGLSestimatorshouldusey=21+.Ifinstead,asin Gryparisetal. ( 2009 ),asmoothedXisderivedusingamixedmodelrepresentationofpenalizedregressionsplinesthatisnotequivalenttokriging,thentherelationshipbetweentheunobservedcovariateandthesmoothedXwillnotfollowaBerksonerrormodel.BecausethekrigingmeanistheBLUP,anyotherestimateofE(X(su)jX(so))(e.g.,oneobtainedfromsmoothing)willnecessarilybebiasedand,whenusedtoaligndisparatedatasets,willintroduceclassicalandBerksonerror. 49

PAGE 50

Toobtainanestimatedgeneralizedleastsquaresestimator,weiteratebetweenestimatingthecovarianceparametersandestimating.Because1ispartoftheinducedmodel'scovariancestructure,estimatingtheparametersischallenging.Inpreviouswork,thecovarianceparameterswereestimatedusingtheOLSresiduals,R=Y)]TJ /F16 11.955 Tf 13.09 0 Td[(X(X0X))]TJ /F5 7.97 Tf 6.58 0 Td[(1X0Y,butassumingtheerrorstructuredenedby( Madsenetal. ( 2008 ); Madsen ( 2004 )).TheWLSestimatoristhencalculatedbyusing^=^R)]TJ /F6 11.955 Tf 13.56 2.65 Td[(^21,OLS.However,theinducedcovariancey=21+isthecorrectcovarianceassociatedwithY)]TJ /F16 11.955 Tf 12.41 0 Td[(X.Ify,thentheestimatesbasedonOLSresidualsmayperformwell.However,thisisafairlyrestrictiveassumptionthatmaynotbesatisedinmanypracticalapplications.Morerecently, Szpiroetal. ( 2011 )proposedestimatingthecovarianceparametersusingaleastsquaresttotheempiricalbinnedsemivariogram.Theempiricalbinnedvariogramiscommonlyusedtoapproximatethetheoreticalvariogrammodel, 2(s1,s2)=Var(Y(s1))]TJ /F4 11.955 Tf 11.96 0 Td[(Y(s2)),(3)asafunctionofdistance.Toderivethevariogramasafunctionofthedistancebetweentwopointsconsider 2(s+h,s)=Var(Y(s+h))]TJ /F4 11.955 Tf 11.96 0 Td[(Y(s)).(3)Forexample,suppose=[2yy]areparametersthatdeneanexponentialcovariancestructurewithscale2yandrangey.Moreover,assumex=[2xx]areparametersthatdeneanexponentialcovariancestructurewithscale2xandrangex.RecalltheinducedmodelY(su)=01n1+1k(su)+,N(0,y=21+).ThetheoreticalvariogramVar(Y(s+h))]TJ /F4 11.955 Tf 12.19 0 Td[(Y(s))foranarbitrarylocationsanddistancehiscalculatedasfollows.Lets+h,sodenotethe1nomatrixthatrepresentsthe 50

PAGE 51

covariancebetweenX(s+h)andtheobservedX(so).Lets,sobedenedsimilarly. Var(Y(s+h))]TJ /F21 10.909 Tf 10.91 0 Td[(Y(s))=Var(Y(s+h))+Var(Y(s)))]TJ /F20 10.909 Tf 10.91 0 Td[(2Cov(Y(s+h),Y(s)) (3) =21[2x)]TJ /F25 10.909 Tf 10.9 0 Td[(s+h,so)]TJ /F5 7.97 Tf 6.58 0 Td[(1oo0s+h,so]+2y+21[2x)]TJ /F25 10.909 Tf 10.91 0 Td[(s,so)]TJ /F5 7.97 Tf 6.58 0 Td[(1oo0s,so]+2y)]TJ /F20 10.909 Tf 8.48 0 Td[(221h2xe)]TJ /F19 5.978 Tf 10.28 3.25 Td[(h x)]TJ /F25 10.909 Tf 10.91 0 Td[(s+h,so)]TJ /F5 7.97 Tf 6.59 0 Td[(1oos,soi+2ye)]TJ /F19 5.978 Tf 10.34 3.25 Td[(h y (3) =2212x[1)]TJ /F21 10.909 Tf 10.91 0 Td[(e)]TJ /F19 5.978 Tf 10.28 3.25 Td[(h x]+22y[1)]TJ /F21 10.909 Tf 10.91 0 Td[(e)]TJ /F19 5.978 Tf 10.33 3.26 Td[(h y])]TJ /F24 10.909 Tf 8.48 0 Td[(21s+h,so)]TJ /F5 7.97 Tf 6.59 0 Td[(1oo0s+h,so)]TJ /F24 10.909 Tf 8.48 0 Td[(21s,so)]TJ /F5 7.97 Tf 6.58 0 Td[(1oo0s,so+221s+h,so)]TJ /F5 7.97 Tf 6.59 0 Td[(1oo0s,so (3) Becausethenalthreetermsin 3 dependonthelocationsofthetwopoints,thevarianceofthedifferencebetweentwolocationsseparatedbyadistancehdependsonhandthelocationofthetwopointsinrelationtothelocationsso.Asaresult,thetheoreticalvariogramdoesnotdependonhaloneandthus,empiricalbinnedvariogramsthatassumestationarityandisotropydonotapproximatethetheoreticalvariogramderivedhere.However,undercertaincombinationsof1,xand2x, 3 mayapproximatelyequal22y[1)]TJ /F4 11.955 Tf 12.45 0 Td[(e)]TJ /F19 5.978 Tf 10.34 3.26 Td[(h y].Thus,insomecases,estimatesbasedonanempiricalbinnedsemivariogrammayperformwell.Theprimarycontributionofthisworkisacarefulconsiderationofhowtoestimateallcovarianceparameterssothatuncertaintyisproperlyquantiedwhendataarespatiallymisaligned.Thisisacrucialsteptoobtaininganestimatedgeneralizedleastsquares(EGLS)estimator.AlthoughtheGLSestimatorinSection 3.3 isbothintuitiveandpractical,detailsregardingitsimplementationandafullevaluationofitsperformanceinthecontextofspatialmisalignmenthaveyettobeprovided.Inthenexttwosections,twoEGLSmethodsforestimatingwhilefullyaccountingforthecovariancestructureoftheinducedmodelaredeveloped. 51

PAGE 52

3.3.1.1Method-of-momentsHere,amethod-of-momentsestimatorofyandthesubsequentestimatedgeneralizedleastsquaresestimatorarederivedwhen=2Inn.First,theOLSestimatorof, ^OLS=(X0X))]TJ /F5 7.97 Tf 6.59 0 Td[(1X0Y(su)(3)isunbiasedfor.Second,theexpression Q=Y0(I)]TJ /F10 11.955 Tf 11.96 0 Td[(P)Y)]TJ /F4 11.955 Tf 11.96 0 Td[(trace((I)]TJ /F10 11.955 Tf 11.96 0 Td[(P)(21)) n)]TJ /F4 11.955 Tf 11.96 0 Td[(p,(3)whereP=X(X0X))]TJ /F5 7.97 Tf 6.58 0 Td[(1X0,=uu)]TJ /F10 11.955 Tf 12.02 0 Td[(uo)]TJ /F5 7.97 Tf 6.58 0 Td[(1ooou,andpisthedimensionof,isunbiasedfor2.Consider ^Q=Y0(I)]TJ /F10 11.955 Tf 11.96 0 Td[(P)Y)]TJ /F4 11.955 Tf 11.95 0 Td[(trace((I)]TJ /F10 11.955 Tf 11.96 0 Td[(P)(^21,OLS)) n)]TJ /F4 11.955 Tf 11.96 0 Td[(p(3)asanestimateofQ.BasedontheOLSestimateof1and^Q,consider ^y=^21,OLS+^QInn(3)asanestimateofy.Finally,anEGLSestimatoris ^EGLS=(X0^)]TJ /F5 7.97 Tf 6.59 0 Td[(1yX))]TJ /F5 7.97 Tf 6.59 0 Td[(1X0^)]TJ /F5 7.97 Tf 6.59 0 Td[(1yY(su),(3)andanestimateofthevarianceis ^Var(^EGLS)=(X0^)]TJ /F5 7.97 Tf 6.59 0 Td[(1yX))]TJ /F5 7.97 Tf 6.59 0 Td[(1.(3)Becauseyisafunctionof1,theinitialestimatesareusedinaniterativeproceduretoobtainthenalparameterestimates: 1. Update^Q. 2. Update 3 with^Qfrom2and^1,EGLS. 3. Recalculate^EGLSanditsvariance,Var(^EGLS). 52

PAGE 53

4. Repeatsteps1-3untilchangesintheestimatearesufcientlysmall. 3.3.1.2RestrictedmaximumlikelihoodbasedmethodsThemethodofmomentsapproachprovidesaframeworkforestimatinginthespecialcasewhen=2.IfthecovariatesandXaccountforallthespatialvariationinY,thentheindependenceassumptionisvalid.Otherwise,Ywouldrequireamorecomplicatederrorstructure.Here,themoregeneralcasewhereisavectorofmparametersisconsidered.Inmixedmodels,thevariancecomponentsofthemodelcanbeestimatedbasedonmaximizingthelikelihoodofameanzeroerrorcontrastinsteadofestimatingthemeanandcovarianceparametersjointlyusingthefulllikelihood.Thatis,thelikelihoodofatransformationofthedatadenedasK0Y(su)ismaximized,whereKischosensuchthatE(K0Y(su))=0.FortheBerksonerrormodelunderconsideration,aniterativealgorithmisrequiredbecause1isapartofboththemeanandthevariance.Thatis,thelikelihoodofameanzerocontrastisderivedandmaximizedusingworkingestimatesof1asiftheywerethetruevaluesof1.Itisworthmentioningthatthetermrestrictedlikelihoodisoftenreservedforanerrorcontrastthatiscompletelyindependentofthemeanparameters.Becausethemeanparameter1isstillapartoftheerrorcontrast'slikelihood,thetermrestrictedlikelihoodisnotused.Instead,weusethetermREML-like.ThefollowingalgorithmisusedtoobtaintheREML-likeestimatesof: 1. EstimateusingOLS ^OLS=(X0X))]TJ /F5 7.97 Tf 6.59 0 Td[(1X0Y(su).(3) 2. EstimatebymaximizingtheloglikelihoodusingtheOLSestimateof1, ^(0)=argmaxl(0)r(3) 53

PAGE 54

wherethelikelihoodl(0)ris )]TJ /F6 11.955 Tf 11.95 0 Td[(2l(0)r=c+log(j+^21,OLSj)+log(jX0(+^21,OLS))]TJ /F5 7.97 Tf 6.59 -.01 Td[(1Xj)+r0(+^21,OLS))]TJ /F5 7.97 Tf 6.58 -.01 Td[(1r (3) wherecisaconstantthatdoesnotdependonor,and r=Y)]TJ /F16 11.955 Tf 11.96 0 Td[(X(X0(+^21,OLS))]TJ /F5 7.97 Tf 6.59 0 Td[(1X))]TJ /F5 7.97 Tf 6.59 0 Td[(1X0(+^21,OLS))]TJ /F5 7.97 Tf 6.58 0 Td[(1Y(su).(3) 3. Basedon^(0),calculate^(0)y=(^(0))+^21,OLSandtheinitialEGLSestimateof, ^(0)EGLS=(X0(^(0)y))]TJ /F5 7.97 Tf 6.58 0 Td[(1X))]TJ /F5 7.97 Tf 6.59 0 Td[(1X0(^(0)y))]TJ /F5 7.97 Tf 6.58 0 Td[(1Y(su)(3) 4. Calculate^(k)EGLSasfollows ^(k)]TJ /F5 7.97 Tf 6.58 0 Td[(1)=argmaxl(k)]TJ /F5 7.97 Tf 6.59 0 Td[(1)rwhere )]TJ /F6 11.955 Tf 11.96 0 Td[(2l(k)]TJ /F5 7.97 Tf 6.59 0 Td[(1)r=c+log(j+(^(k)]TJ /F5 7.97 Tf 6.58 0 Td[(1)1,EGLS)2j)+log(jX0(+(^(k)]TJ /F5 7.97 Tf 6.59 0 Td[(1)1,EGLS)2))]TJ /F5 7.97 Tf 6.59 0 Td[(1Xj)+r0(+(^(k)]TJ /F5 7.97 Tf 6.59 0 Td[(1)1,EGLS)2))]TJ /F5 7.97 Tf 6.59 0 Td[(1r (3) and r=Y)]TJ /F16 11.955 Tf 11.95 0 Td[(X(X0(+(^(k)]TJ /F5 7.97 Tf 6.59 0 Td[(1)1,EGLS)2))]TJ /F5 7.97 Tf 6.59 0 Td[(1X))]TJ /F5 7.97 Tf 6.58 0 Td[(1X0(+(^(k)]TJ /F5 7.97 Tf 6.59 0 Td[(1)1,EGLS)2))]TJ /F5 7.97 Tf 6.58 0 Td[(1Y(su)(3)Basedon^(k)]TJ /F5 7.97 Tf 6.59 0 Td[(1),calculate^(k)]TJ /F5 7.97 Tf 6.58 0 Td[(1)y=(^(k)]TJ /F5 7.97 Tf 6.59 0 Td[(1))+(^(k)]TJ /F5 7.97 Tf 6.58 0 Td[(1)1,EGLS)2andthesubsequentnewEGLSestimateof ^(k)EGLS=(X0(^(k)]TJ /F5 7.97 Tf 6.58 0 Td[(1)y))]TJ /F5 7.97 Tf 6.58 0 Td[(1X))]TJ /F5 7.97 Tf 6.59 0 Td[(1X0(^(k)]TJ /F5 7.97 Tf 6.59 0 Td[(1)y))]TJ /F5 7.97 Tf 6.59 0 Td[(1Y(su)(3) 5. Repeatstep 4 untilchangesin^EGLSand^aresufcientlysmall.Then,anestimateofthevarianceofthenalestimator,^EGLS,is ^Var(^EGLS)=(X0^)]TJ /F5 7.97 Tf 6.59 0 Td[(1yX))]TJ /F5 7.97 Tf 6.59 0 Td[(1.(3) 54

PAGE 55

3.3.2EstimatingTheBerksonParametersRecall,theBerksonparameters=f,xgwereassumedknowntoderivetheMOMbasedEGLSestimatoranderrorcontrastmaximumlikelihoodbasedEGLSestimatorinSections 3.3.1.1 and 3.3.1.2 ,respectively.Forsomecasesinpractice,Berksonparametersareestimatedusingvalidationdata.Here,theobservedX(so)canbeconsideredvalidationdataandandxcanbeestimatedusingrestrictedmaximumlikelihoodandEGLS.Thatis,^xisfoundbymaximizingtherestrictedlikelihoodofxjX(so).Then, ^EGLS=(C(so)0^ooC(so)))]TJ /F5 7.97 Tf 6.59 0 Td[(1C(so)0^ooX(so),(3)where^oo=oo(^x).If^=f^EGLS,^xgisusedinthelikelihoodofY(su)jX(so),Z,,,asiftheyareknown,i.e.,Y(su)jX(so),Z,,,^,thenthemethodsdescribedinSections 3.3.1.1 and 3.3.1.2 canbeused.Inthecontextofthemeasurementerrorliterature,theestimatorsarebasedonpseudo-likelihoodapproaches.Alimitationofpseudo-likelihoodapproachesisestimatingintroducesacomponentofclassicalmeasurementerrorthatcancausebiasinparameterestimatesandstandarderrorstobetoosmall( Buonaccorsi ( 2009 )).Bootstrapmethodsarealternativestothepseudo-likelihoodapproachesthatcanaccountfortheBerksonandclassicalmeasurementerror( Szpiroetal. ( 2011 )).However,thesemethodsarecomputationallyintensiveandifthenumberofobservationsinX(so)islarge,thentheeffectofclassicalmeasurementerrorisnotmuchofanissue( Buonaccorsi ( 2009 )).Nonetheless,itisusefultoaccountfortheclassicalmeasurementerrorinducedfromestimating.Recall, 264X(su)X(so)375j,xN0B@264C(su)C(so)375,x(x)=264uu(x)uo(x)ou(x)oo(x)3751CA,(3) 55

PAGE 56

and X(su)jX(so),,xN(X(su)=C(su)+uo)]TJ /F5 7.97 Tf 6.59 0 Td[(1oo(X(so))]TJ /F10 11.955 Tf 11.95 0 Td[(C(so)),=uu)]TJ /F10 11.955 Tf 11.96 0 Td[(uo)]TJ /F5 7.97 Tf 6.58 0 Td[(1ooou). (3) Toremoveconditioningon,considerauniformpriorfor.Inthiscase,theposteriordistributionjX(so),xis jX(so),xN(~GLS=(C(so)0)]TJ /F5 7.97 Tf 6.58 0 Td[(1ooC(so)))]TJ /F5 7.97 Tf 6.59 0 Td[(1C(so)0)]TJ /F5 7.97 Tf 6.59 0 Td[(1ooX(so),(C(so)0)]TJ /F5 7.97 Tf 6.58 0 Td[(1ooC(so)))]TJ /F5 7.97 Tf 6.59 0 Td[(1). (3) Theconditionaldensityis (X(su)jX(so),x)=Z(X(su)jX(so),x,)(jX(so),x)d, (3) wherethenotation(xjy)representsthedensityofxgiveny.Basedontheintegral 3 X(su)jX(so),xN(~X(su)=C(su)~GLS+uo)]TJ /F5 7.97 Tf 6.59 0 Td[(1oo(X(so))]TJ /F10 11.955 Tf 11.96 0 Td[(C(so)~GLS),~=uu)]TJ /F10 11.955 Tf 11.96 0 Td[(uo)]TJ /F5 7.97 Tf 6.58 0 Td[(1ooou+(C(su))]TJ /F10 11.955 Tf 11.96 0 Td[(uo)]TJ /F5 7.97 Tf 6.59 0 Td[(1ooC(so))(C(so)0ooC(so)))]TJ /F5 7.97 Tf 6.59 0 Td[(1(C(su))]TJ /F10 11.955 Tf 11.96 0 Td[(uo)]TJ /F5 7.97 Tf 6.59 0 Td[(1ooC(so))0). (3) Now,theBerksonparametersarexandthemethodsinSections 3.3.1.1 and 3.3.1.2 canbereformulatedbyreplacingX(su)andwith~X(su)and~,respectively.Thenthepseudo-likelihoodapproachesinSection 3.3.2 wouldonlyfailtoincorporatetheclassicalmeasurementerrorassociatedwithestimatingx. 56

PAGE 57

3.4Applications 3.4.1EPAEnvironmentalMonitoringandAssessmentProgramTheU.S.EnvironmentalProtectionAgency(EPA)conductedtheEnvironmentalMonitoringandAssessmentProgramfrom1990-2006.Theprogram'smaingoalwastocollectdatatofacilitatetheadvancementoftheunderstandingofecologicalrisk.Morespecically,onearmoftheprogramcollecteddataonthehealthofaquaticecologicalresourcestoestimatecurrentandfuturerisks( U.S.EnvironmentalProtectionAgency ( 1999a ); U.S.EnvironmentalProtectionAgency ( 1999b )).Asin Szpiroetal. ( 2011 )and Madsenetal. ( 2008 ),weconsidertherelationshipbetweenchloridelevelsinstreamsandthewatershedforestationlevel.Followingpreviouswork,weusethenaturallogarithmofthechlorideconcentrationsastheresponsevariable,Y,andatransformationofthewatershedforestationlevel,logit((.98%Forested+1)=100),asthecovariate,X.Thedata,collectedbetween1993and1996,arefromtheMid-AtlanticHighlandsregionoftheeasternUnitedStates.Whenmultiplemeasurementsfromdifferenttimeswerecollectedatthesamelocation,weusetheearliesttime( Szpiroetal. ( 2011 )).Bothresponseandcovariateareobservedat422locations.Atanadditional157locations,onlytheresponseisobserved(Figure 3-1 ).Forthealigneddata,thatis,the422observationswhereboththeresponseandcovariateareobserved,weusedrestrictedmaximumlikelihoodtotasimplelinearregressionmodelrelatingtheresponsetothecovariateassuminganexponentialcovariancestructuredenedbythenugget,scale,andrange.Basedonthist,wefoundahighlysignicantnegativeassociation,^1=)]TJ /F6 11.955 Tf 9.3 0 Td[(0.32,SE=0.019.Forthecaseofthemisaligneddata,werstusedthe422observationsofthecovariateandrestrictedmaximumlikelihoodtotalinearmodelwithaconstantmeanandanerrorstructuredenedbyanexponentialcovariancestructure.Anexponential 57

PAGE 58

covariancestructurewaschosenbasedonvisuallyassessinganempiricalbinnedsemivariogramofthealigneddata.Moreover,theexponentialcovariancestructurewasalsousedin Szpiroetal. ( 2011 ).Theestimatedmeanparameter^was1.41andtheestimatedcovarianceparameters^xwerenugget=0.768,scale=4.378,range=23.841.ThekrigingdistributionofXatthe157locationswhereYisobservedwasusedtorecalculatetheregressioncoefcient1using 1. TheerrorcontrastmaximumlikelihoodmethoddescribedinSection 3.3 2. OrdinaryLeastSquares(OLS)assuminganindependenterrorstructure 3. REMLandEstimatedGeneralizedLeastSquaresassuminganexponentialcovariancestructure 4. RobustregressionusingR'srlmfunction 5. Theparametricbootstrap(PB)( Szpiroetal. ( 2011 ))with1000bootstrapsamples 6. Thepartialparametricbootstrap(PPB)( Szpiroetal. ( 2011 ))1000bootstrapsamples 7. AWLSmethod( Madsenetal. ( 2008 ))Noteforthebootstrapmethods,thecovarianceparameterswereestimatedusingtheREMLmethoddescribedinSection 3.3 .TheresultsareprovidedinTable 3-1 .TheWLSestimatoryieldedasubstantiallydifferentregressionestimate.AsshowninSection 3.3.1 ,theestimateofusedastheweightingmatrixdependsonanestimateof1.Afterfurtherinspection,wediscovered,inthiscase,theWLSestimatewassensitivetotheOLSestimateof1usedtoestimate^.Asanexample,iftheinitialestimate1wasincreasedby0.001and0.01,theWLSestimateswouldbe-0.42and-0.17,respectively,withstandarderrorsof0.052and0.045.TheregressioncoefcientestimatedusingOLS,PB,PPB,androbustregressionaresimilar,)]TJ /F6 11.955 Tf 9.3 0 Td[(0.37,)]TJ /F6 11.955 Tf 9.29 0 Td[(0.37,)]TJ /F6 11.955 Tf 9.3 0 Td[(0.37and)]TJ /F6 11.955 Tf 9.3 0 Td[(0.40respectively.Therobustregression,PB, 58

PAGE 59

andPPBestimatesaresimilartotheestimateobtainedby Szpiroetal. ( 2011 ).However,afteraccountingfortheinducedcovarianceusingourmethodorbyttingaspatiallinearmodel,theestimatedregressioncoefcientdecreasessubstantiallyto)]TJ /F6 11.955 Tf 9.3 0 Td[(0.25and)]TJ /F6 11.955 Tf 9.3 0 Td[(0.27forthespatialGLSestimatorandthenewREML-likemethodweproposed,respectively.TheestimatedstandarderrorsbasedonourmethodandthespatiallinearmodelareslightlygreaterthantheOLSstandarderrorthatassumesindependenterrorsandthestandarderrorobtainedfromrobustregression.ThestandarderrorsobtainedfromusingthePBandPPBarelargerthananyoftheothermethodswhichisconsistentwithsimulationstudies( Lopianoetal. ( 2011 ); Szpiroetal. ( 2011 )).Inthisexample,accountingfortheinducedcovariancestructurehasasubstantialeffectontheparameterestimate(seeFigure 3-2 ).ThePBandPPBpointestimatesdonotidentifythischangeastheyuseOLStoestimate.Notethattheestimatedslopebasedonthemisaligneddataandtheestimatedslopeusingthealigneddataarewithinastandarderrorofeachother.Althoughthepoint-to-pointmisalignmentproblemiscommoninmanyecologicalstudies,thepoint-to-arealmisalignmentproblemismorecommoninpublichealthassociationstudies.Forcondentialityreasons,publichealthdataareusuallyaggregatedtothearealunitlevel(e.g.,county)whileenvironmentalriskfactorsareobservedatpoints.Inthenextexample,weextendourmethodologytothecasewhentheresponseisobservedatthearealunitlevelandthecovariateisobservedatthepointlevel. 3.4.2CDCEnvironmentalPublicHealthTrackingProgramOnegoaloftheEnvironmentalPublicHealthTracking(EPHT)Programstudyistoassesstheassociationbetweenairquality(here,ozone)andmyocardialinfarctions(MI).OzonemeasurementsarerecordedfromanetworkofairmonitorsplacedthroughoutthestatebytheFloridaDepartmentofEnvironmentalProtection(Figure 3-3 ).DuringthestudyperiodofAugust2005,ozonedatawereavailablefrom56monitors.Themaximumofthedailymaximum8-houraverageozonevaluesduringamonthisused 59

PAGE 60

asthemonthlydatavalueforaparticularmonitorandisreferredtoasthemonthlymaximumozone.Florida'sDepartmentofHealth(FDOH)hasadata-sharingagreement(DSA)withFlorida'sAgencyforHealthCareAdministration(AHCA),whichprovidedaccesstotwodatasources:condentialhospitalizationrecordsandemergencyroomrecords.ThehospitalizationrecordscontainalladmissionstoFlorida'spublicandprivatehospitalswhereeithertheprimaryorsecondarycauseofadmissionwasMI(ICD-9codes410.0.0)( NationalCenterforHealthStatistics ( 1998 )).Ifapersonpresentstotheemergencyroomandissubsequentlyadmittedtothehospital,thatperson'srecordisincludedinthehospitalizationrecords,butnottheemergencyroomrecords.Non-Floridaresidentswereexcludedfromtheanalysis.ConsistentwithourDSA,AHCAprovidedboththezipcodeandcountyofresidenceforeachpatient'srecord.Selectedpatientdemographicinformationisalsorecorded,includingsex,age,andrace/ethnicity.Foranalysis,thedatamustbelinkedonthesamegeographicalunits,which,forcondentialityreasons,hasbeenthecountylevel( Youngetal. ( 2009a ); Youngetal. ( 2009b )).TheMIdatawereindirectlystandardizedbyage(age45,45,55,and>65years),race/ethnicity(black,white,orother),andsex(femaleormale).Theinformationregardingage,race/ethnicity,andsexfortheMIcaseswasobtainedfromthehospitalrecords.Florida'spopulationwasusedasthecomparisonstandardtocalculatethenumberofexpectedMIcases.ThisgaveanMIstandardizedeventratio(MISER),denedastheratioofthenumberofobservedMIcasestothatexpectedamongtheFloridapopulation,foreachcounty( Woodward ( 2004 )).LetY(B)bethenaturallogarithmoftheMISERinthe67countiesinFlorida,wherethecountiesarerepresentedbyB.Becauseozone,X,isobservedatthepointlevel,theblockkrigingdistributionofozoneatthearealunitleveliscalculatedbasedonassuminganexponentialcovariancestructuredenedbytheparametersxwithanexposuresurface 60

PAGE 61

trenddenedby, X(B)jX(so),,xN(X(B)=C(B)+B,so)]TJ /F5 7.97 Tf 6.59 0 Td[(1so(X(so))]TJ /F10 11.955 Tf 11.96 0 Td[(C(so)),=B,B)]TJ /F10 11.955 Tf 11.96 0 Td[(B,so)]TJ /F5 7.97 Tf 6.59 0 Td[(1soso,B) (3) wheresocorrespondstothemonitorlocations,C(B)isamatrixwherethei,jthentryistheaverageofthejthexposuresurfacecovariatewithinblocki,i=1,...,n,B,Bisthecovariancebetweentheblocks,B,soisthepoint-to-blockcovarianceand,soisthecovariancematrixamongthepoints.Thatis, X(B)jX(so),,x=X(B)+v,vN(0,)(3)NoticetherelationshiptakesthesameBerksonformasinthepoint-to-pointmisalignmentproblem.Asaresult,thesamemethodscanbeusedtoestimatetheparametersinthemodel log(Y(B))=01671+1X(B)+Z(B)[2,...,6]0+,N(0,2I6767)X(B)jX(so),,x=X(B)+v,vN(0,) (3) whereZ(B)isamatrixofcovariateswhereeachrowcorrespondstothecovariatesforaparticularcounty.Blockkrigingestimateswerederivedassumingaconstantmeanwithanexponentialcovariancestructure( Lopianoetal. ( 2011 )).Theestimatedblockkrigingmeanswereusedinthemodelin( 3 ).Theestimatedregressioncoefcientassociatedwithozonewas0.025withastandarderrorof0.018.Inthisstudy,ozoneisnotasignicantpredictorofthenaturallogarithmoftheMISER.TheestimatedregressioncoefcientandstandarderrorwerethesameforboththeREML-likeandMOMbasedestimates.TheestimatecalculatedusingOLSwasapproximatelyequaltotheadjustedestimateswithaslightlylargerstandarderror(Table 3-2 ).Afterfurtherinspection,weseetheestimateof2derivedwithoutaccountingforthekrigingvariancetermislargerthantheestimatesbasedonREML-likeandMOM 61

PAGE 62

(Table 3-3 ).Theincreaseintheestimateof2isduetotheadditionalvariationinducedfromthekrigingvariance.Thatis,whennotaccountingfortheadditionalvariationfromthekrigingvariance,theestimateof2istheMSE.TheexpectedvalueoftheMSEisgreaterthan2(seeSection 3.3.1.1 ). 3.5Simulation 3.5.1Point-to-PointMisalignmentAsimulationstudywasusedtoassesstheperformanceoftheproposedmethodologiesinthepoint-to-pointmisalignmentsetting.ThedesignofthestudyisbasedontheEPAapplicationandthesimulationstudyof Szpiroetal. ( 2011 ).TheresponsevariableYwasgeneratedfromthemodelin 3 with0=5,1=)]TJ /F6 11.955 Tf 9.3 0 Td[(0.321andnoadditionalcovariatesZ.WeconsideredtwoerrorstructuresforY,independentandcorrelated.Fortheindependentcase,N(0,0.858I).Forthecorrelatedcase,wasdenedbyanexponentialcovariancestructurewithnugget,scaleandrangedenedastherespectiveelementsof=[0.496,0.8167,198].ThejointdistributionofX(su)andX(so)isspeciedasin 3 withaconstantmean=1.408andanexponentialcovariancestructurewiththenugget,scaleandrangedenedastherespectiveelementsofthecovarianceparametervectorx=[0.768,4.378,23.841].Weassumethatthecovariateisobservedat422randomlychosenlocationsina300500domainwhiletheresponseisobservedat157randomlychosenlocationsinthesamedomain(Figure 3-4 ).Theresultsfrom3000MonteCarlosimulationsarefoundinTable 3-4 .Theaverage,^1,andvariance,^2^1,of^1fromthesimulateddatasets,andtheaverageestimatedvarianceof^1,s2^1,areprovided.If^1isanunbiasedestimatorof1then^1shouldbewithinerrorofthetruevalueof1=)]TJ /F6 11.955 Tf 9.29 0 Td[(0.321.Ifthevariances2^1of^1isestimatedunbiasedly,s2^1,theaverageoftheestimatedvariancesof^1and^2^1shouldbeapproximatelyequal. 62

PAGE 63

Forthecaseofindependentoutcomes(Table 3-4 ,Rows1-6),theREML-likeandMOMestimatorsperformedsimilarly.WhentheBerksonparameterswereassumedknown,thepointestimateof1wasnearlyunbiasedandthevariabilitywasslightlyover-estimatedonaverage.Asaresult,coveragefornominal95%condenceintervalswasapproximately96%forbothREML-likeandMOM.WhentheBerksonparameterswereestimated,nearlyunbiasedestimatesof1wereobtained.Again,thevariabilitywasslightlyover-estimatedandasaresult,coveragefor95%condenceintervalswas95.4%and95.2%forREML-likeandMOMrespectively.Here,whentheeffectduetomisalignmentisignored(Table 3-4 ,Rows5-6,MethodKrigeandRegress(KR)),i.e.,thekrigingmeanisusedasifitisthetruecovariate,nearlyunbiasedestimatesoftheparameterandtheuncertaintyareobtained.ForknownandunknownBerksonparameters,coveragefor95%condenceintervalsarenear95.5%and95.3%,respectively.Here,anindependenterrorstructureisassumedwhenusingtheKRmethod.Thus,insomecases,theeffectduetomisalignmentcanbeignored.However,futureresearchisneededtodeterminewhenthisistrue.Forthecaseofcorrelatedoutcomes(Table 3-4 ,Rows7-10),whentheBerksonparameterswereknown,1wasestimatedunbiasedlyandthevariabilitywasslightlyunderestimated.Asaresult,coveragefor95%condenceintervalswas94.6%.WhentheBerksonparameterswereestimated,theparameterestimatewasnearlyunbiasedandthevariabilitywasslightlyunder-estimated.Asaresult,coveragefornominal95%condenceintervalswasapproximately94%.Unlikeinthecaseofindependentresponses,theKRmethodunderestimatesthevariabilityandcoveragefor95%condenceintervalsisnear80%.Here,anexponentialcovariancestructureisassumedwhenusingtheKRmethod.TheperformanceofthenaivemethodsherefurtherillustratestheneedformoreresearchtodeterminewhenBerksonerrorsubstantiallyaffectstheperformanceofmethodsthatignoretheBerksonerror. 63

PAGE 64

3.5.2Point-to-ArealMisalignmentAsecondsimulationstudywasusedtoassesstheperformanceoftheproposedmethodologiesinthepoint-to-arealmisalignmentsetting.AsshowninSection 3.4.2 ,themethodsforthepoint-to-pointmisalignmentproblemarealsovalidforproblemswhereblockkrigingisusedtopredictacovariateobservedatthepointlevelatthearealunitlevel.ThedesignofthestudyisbasedontheEPHTapplicationandthesimulationin Lopianoetal. ( 2011 ).Twodatageneratingmechanismsfortheresponsevariableareconsidered.Forbothresponsedatageneratingmechanisms,thejointdistributionofXatthelocationsinFigure 3-5 ismultivariatenormalwithaconstantmean=40andanexponentialcovariancestructurewiththescaleandrangedenedastherespectiveelementsofthecovarianceparametervectorx=[70,200].Intherstdatageneratingsetting,theobservationsaregeneratedatthepointlevel(Figure 3-5 )accordingtoalinearmodelwith0=)]TJ /F6 11.955 Tf 9.3 0 Td[(0.8,1=0.2,noadditionalcovariatesZandN(0,42I).Then,thearealunitresponseistakentobetheaverageofthepointlevelresponses.Intheseconddatageneratingsetting,thearealunitresponseisgeneratedbyrstaveragingthecovariatevalueswithinthearealunitandthengeneratingtheresponseatthearealunitlevelaccordingtolinearmodelwith0=)]TJ /F6 11.955 Tf 9.3 0 Td[(0.8,1=0.2andnoadditionalcovariatesZandN(0,42I).Weassumethatthecovariateisobservedatthe56monitorlocationsshowninFigure 3-3 whiletheresponseisobservedatthe67countiesinFlorida.Foreachmethod,blockkrigingmeansareusedtoalignthedatasetsatthecountylevelandthenregressioniscarriedoutbytakingintoaccounttheadditionalvariationinducedbytheblockkrigingvariance.Theresultsfrom3000MonteCarlosimulationsarefoundinTable 3-5 .Theaverage,^1,andvariance,^2^1,of^1fromthesimulateddatasets,andtheaverageestimatedvarianceof^1,s2^1,areprovided.Asinthelastexample,if^1isanunbiasedestimatorof1then^1shouldbewithinerrorofthetruevalueof1=0.2.Ifthevariances2^1of^1is 64

PAGE 65

estimatedunbiasedly,s2^1,theaverageoftheestimatedvariancesof^1and^2^1shouldbeapproximatelyequal.Forthecasewhentheresponseisgeneratedatthepointlevel,theMOMandREML-likeestimatorsof1areunbiasedandprovidenearlyunbiasedestimatesofuncertaintywhentheBerksonparametersareknown.Asaresult,coveragefor95%condenceintervalsis95.3%forbothcases.WhentheBerksonparametersareestimated,theMOMandREML-likeestimatorsareslightlybiasedbuttheestimatesofuncertaintyareessentiallyunbiased.Coveragefor95%condenceintervalsare94.7%and94.9%fortheMOMandREML-likeestimators,respectively.Finally,thekrigeandregress(KR)methodthatassumesindependenterrorstructuresprovidesanunbiasedestimateoftheparameterandanessentiallyunbiasedestimateoftheuncertaintywhentheBerksonparametersareknown.Coveragefor95%condenceintervalsis95%.However,theuncertaintyintheparameterestimateisalmostthreetimesaslargeastheothermethods.ForthecasewhentheBerksonparametersareunknown,theKRestimatorisslightlybiasedupwardandtheestimateofuncertaintyisbiaseddownward.Asaresultcoveragefor95%condenceintervalsisnear95%.Forthecasewhentheresponseisgeneratedatthearealunitlevel,theMOMandREML-likeestimatorsareunbiasedandprovideunbiasedestimatesofuncertainty.Coveragefor95%condenceintervalsare94.9%and95.1%fortheMOMandREML-likeestimators,respectively.WhentheBerksonparametersareestimated,theMOMandREML-likeestimatorsareslightlybiasedbuttheestimatesofuncertaintyareessentiallyunbiased.Coveragefor95%condenceintervalsare94.9%and94.9%fortheMOMandREML-likeestimators,respectively.Here,theKRmethodyieldedthesameresultsasinthecasewheretheresponseisgeneratedatthepointlevel.However,inthiscase,theuncertaintyintheparameterestimateismorethan30timesaslargeastheothermethods. 65

PAGE 66

3.6ConcludingRemarksInthespatialmisalignmentframeworkconsideredhere,acomplicatederrorstructurethatinvolvestheunknownparameterofinterestarisesmakingcovarianceparameterestimationadifculttask.Although Szpiroetal. ( 2011 )and Madsenetal. ( 2008 )consideredcovarianceparameterestimation,theirmethodsdonotfullyaccountforthecovariancestructureoftheinducedmodel.Inpractice,estimatedgeneralizedleastsquaresestimatorscannotbeobtainedwithoutestimatesofthecovarianceparameters.OurmaincontributionisadetailedenumerationofhowtoestimatethecovarianceparametersandcalculateanEGLSestimatoroftheregressionparameterofinterest.ByworkingwiththeinducedmodelassumingtheBerksonparametersareknown,weareabletoobtainREML-likeandMOMestimatesofthecovarianceparametersassociatedwiththeoriginalregressionerror.Then,toaccountfortheBerksonerror,theestimatedcovarianceparametersareusedtocalculateanestimatedgeneralizedleastsquaresestimator.Two-stagemethodsarecommonlyusedinpractice( Buonaccorsi ( 2009 )).Moreover,jointmaximumlikelihoodmethods,fullyBayesianmethods,andbootstrapmethodscanbecomputationallyintensiveandthus,two-stagemethodsmaybepreferred.Inthecaseofusingkrigingtoaligndisparatedatasetsinapoint-to-pointmisalignmentproblem,oursimulationresultsindicatetheproposedtwo-stagemethodologyperformswell.Nonetheless,weconsideredabootstrapproceduretoaccountfortheclassicalmeasurementerror.Becauseofthelargenumberoflocationswherethecovariatewasobserved,theeffectduetoclassicalmeasurementerrorwasnegligible.Thatis,whenestimatingtheBerksonparameters,coveragefor95%condenceintervalswasstillapproximately95%.Therefore,theadditionalcomputationaleffortrequiredtoaccountfortheclassicalmeasurementerrorusingabootstrapwasunnecessary.Althoughnotshownhere,abootstrapmethodtoapproximatethesamplingdistributionoftheEGLSestimator,asopposedtotheOLSestimator,caneasilybederived. 66

PAGE 67

Ourresultsextendtothepoint-to-arealmisalignmentproblemwhenblockkrigingisusedtoalignthevariablesatthearealunitlevel.Oursimulationresultsindicatethemethodologyperformswellandisrobusttothedatageneratingmechanism.Inaddition,asin Lopianoetal. ( 2011 ),theeffectofestimatingtheBerksonparametersappearstobenegligibleinsomepoint-to-arealmisalignmentproblems.Becausemostpublichealthassociationstudiesutilizeaggregatedataforcondentialitypurposes,themethodspresentedhereareusefulforroutineassessmentofenvironmentalriskfactors.However,asmoredataareintegratedfrommultiplesources,morecomplicatederrors-in-variablesproblemswillariseandmustbeaddressedinordertoassesstherelationshipsofseveralmisalignedvariables.Thatis,methodsthataccountforallsourcesoferrorareneededtoaccuratelyassesscomplexrelationshipsandmakeinformedpolicydecisions.Support TheworkofthersttwoauthorswassupportedbytheNationalAeronauticalandSpaceAdministration(NASA)AppliedSciencesPublicHealthProgram'sgrantawardnumberNNX10AO08GandbytheFloridaDepartmentofHealth,DivisionofEnvironmentalHealthandGrant/CooperativeAgreementNumber1U38EH000941-01fromtheCentersforDiseaseControlandPrevention(CDC).TheworkoftherstauthorwaspartiallysupportedbytheNationalScienceFoundationunderGrantNo.0801544intheQuantitativeSpatialEcology,EvolutionandEnvironmentProgramattheUniversityofFlorida.ThendingsandconclusionsinthisreportarethoseoftheauthorsanddonotnecessarilyrepresenttheviewsofNASA,CDCorNSF.WewouldliketoacknowledgediscussionswithDr.NoelCressieintheearlystagesoftheresearch.WewanttothankDr.AdamSzpiroforhelpfuldiscussionsandDr.LisaMadsenforhelpfuldiscussionsoftheWLSmethodandforsharinghercodewithus. 67

PAGE 68

Table3-1. EMAPApplication:Estimatedregressioncoefcientandassociatedstandarderrorsbasedondifferentmethods. Method^1StandardError NewREML-like-0.270.063SpatialGLS-0.250.065OLS-0.370.062Robust-0.400.062PB-0.370.092PPB-0.370.091WLS-1.150.064 Table3-2. EPHTApplication:Estimatedregressioncoefcientsandmeasuresofuncertaintybasedondifferentmethods Method^1s2^1s^1 NewREML-like0.0250.000320.018NewMOM0.0250.000320.018OLS0.0240.000330.018 Table3-3. EPHTApplication:Estimatesofregressionerrorbasedondifferentmethods Method^2 NewREML-like0.125NewMOM0.125OLS0.148 68

PAGE 69

Table3-4. Point-to-PointSimulationResults:Average^1andVariance^2^1oftheEstimatesof1,AverageoftheEstimatedVariances2^1of^1,MeanSquaredError(MSE),andCoveragefortheSimulatedDatasetsintheIndependentOutcomeScenario.Thetruevalueof1is)]TJ /F6 11.955 Tf 9.3 0 Td[(0.321. OutcomeMethodBerksonParameters^1^2^1s2^1MSECoverage IndependentREML-likeKnown-0.3220.003000.003170.0030096.0Estimated-0.3230.003120.003220.0031295.4MOMKnown-0.3220.003000.003140.0030095.9Estimated-0.3230.003120.003210.0031295.2KRKnown-0.3220.003010.003100.0030195.5Estimated-0.3230.003150.003150.0031595.3CorrelatedREML-likeKnown-0.3210.003890.003830.0038994.6Estimated-0.3220.004160.003890.0041694.0KRKnown-0.3190.008170.003490.0081880.0Estimated-0.3190.008440.003550.0084480.1 Table3-5. Point-to-ArealSimulationResults:Average^1andVariance^2^1oftheEstimatesof1,AverageoftheEstimatedVariances2^1of^1,MeanSquaredError(MSE),andCoveragefortheSimulatedDatasetsinbothPoint-to-ArealSettings.Theoutcomewasgeneratedatthepointlevelandthearealunitlevel.Thetruevalueof1is0.2 OutcomeMethodBerksonParameters^1^2^1s2^1MSECoverage PointLevelMOMKnown0.2000.001150.001210.0011595.3Estimated0.2040.001350.001340.0013794.7REML-likeKnown0.2000.001150.001200.0011595.3Estimated0.2040.001350.001330.0013694.9KRKnown0.2000.003510.003360.0035194.9Estimated0.2050.004190.003820.0042294.3ArealLevelMOMKnown0.2000.001100.001100.0011094.9Estimated0.2050.001260.001240.0012694.9REML-likeKnown0.2000.001100.001100.0011095.1Estimated0.2050.001260.001250.0012694.9KRKnown0.2000.04220.04160.042294.6Estimated0.2050.04820.04740.048294.5 69

PAGE 70

Figure3-1. ObservationLocations 70

PAGE 71

Figure3-2. Fittedregressionlinesfordifferentmethods. 71

PAGE 72

Figure3-3. MonitorLocationsandCountiesintheStateofFlorida 72

PAGE 73

Figure3-4. LocationsofPointsforSimulationStudy 73

PAGE 74

Figure3-5. GridLocationsandCountiesintheStateofFlorida 74

PAGE 75

CHAPTER4BERKSONERRORINGENERALIZEDLINEARMODELSBriefSummary GeneralizedlinearmodelsareroutinelyusedtoassesstherelationshipbetweenaresponseYandacovariateX.OneassumptionofthemodelisthatthecovariateXismeasuredwithouterror.Becauseinpracticethisassumptionisoftenviolated,errors-in-variablesmodelsweredevelopedtosimultaneouslymodeltherelationshipbetweenYandXandmodeltheerrorinX.TwomodelscommonlyusedtomodeltheerrorinXaretheclassicalmeasurementerrormodelandtheBerksonerrormodel.Inaspatialsetting,whenthegeographicalunitsassociatedwiththeobservedresponsevaluesandobservedcovariatevaluesaredifferent,i.e.,spatiallymisaligned,ameasurementproblemcanarisewhenaligningthedatasets.Ifkrigingisusedtoalignthedisparatedatasets,thenaBerksonerrormodelarises.Withthespatialmisalignmentproblemasamotivatingexample,weconsideraBerksonerrormodelwiththerelationshipbetweenthetruecovariateandresponsedenedbyageneralizedlinearmodel.WeshowhowmethodsusedtoconductinferenceinlinearmixedmodelsandgeneralizedlinearmixedmodelscanbeadaptedtoconductinferenceingeneralizedlinearmodelswithBerksonerror.Thepropertiesoftheinferenceproceduresareassessedtheoreticallyandviasimulationandillustratedusingtwospatialmisalignmentexamplesrelatingairqualityandhealthoutcomes. 4.1BackgroundInformationErrors-in-variablesmodelsareregressionmodelsthataccountformeasurementerrorinoneormorepredictorvariables.Forexample,supposethecovariateXandtheresponsevariableYarerelatedthroughageneralizedlinearmodelwithonecovariate: YExponentialFamily(mean=,scale=)g()=01n1+1X (4) 75

PAGE 76

wheretheithelementsofYandXareassociatedwiththeithobservationalunitandg(),thelinkfunction,isasmoothmonotonicfunction.ThetwomeasurementerrormodelsmostcommonlyusedtodescribetheerrorinXaretheclassicalmeasurementerrormodelandtheBerksonerrormodel( Buonaccorsi ( 2009 ); Carrolletal. ( 2006 )).Inclassicalmeasurementerrormodels,theobservedvalueWgiventhetruevalueXismodeled: W=X+w,(4)wherewisameanzerorandomvariableindependentofXandY.BerksonerrormodelsmodelthetruevalueXgiventheobservedvalueW: X=W+x,(4)wherexisameanzerorandomvariableindependentofWandY.Forbothmodels,iftheconditionaldistributionofYjX,WisequivalenttothedistributionofYjXthenWissaidtobeasurrogateforX.Ifthemomentsordistributionoftheerrormodelaredenedparametrically,thentheparametersarereferredtoasthemeasurementerrorparameters.InBerksonerrorproblems,Wmaydependontheunknownmeasurementerrorparametersandgenerally,inpractice,theparametersandthepredictorWareestimated.Foramoredetaileddiscussionofbothtypesofmeasurementerrormodelssee Buonaccorsi ( 2009 )and Carrolletal. ( 2006 ).HereweconsidertheBerksonerrormodelsinageneralizedlinearmodel(GLM)framework.Thus,themeasurementerrorframeworkisdenedas YExponentialFamily(mean=,scale=)g()=01n1+1XX=W+x,xN(0,x), (4) 76

PAGE 77

Equivalently,wecanwriteY(su)ExponentialFamily(mean=,scale=)g()=01n1+1W+1xxN(0,x). (4)Measurementerrormodelsplayanimportantroleinmanyapplicationsinvolvingspatialdata.Inparticular,inepidemiologicalstudiesrelatingpublichealthoutcomes,suchasstrokesormyocardialinfarctions,tomeasuresofairquality,suchasozoneorparticulatematterconcentrations,thehealthdataandexposuredataareoftencollectedatdifferentsupports,i.e.,spatiallymisaligned.Forexample,ahealtheventmaybelinkedgeographicallytotheresidenceofthepatientwhileanairqualitymeasureisrecordedatamonitoringsite,leadingtopoint-to-pointmisalignment.Morecommonly,healthoutcomesareaggregatedtothearealunitlevel,suchasthecountyorcensusblock( Gelfandetal. ( 2001 ); Lopianoetal. ( 2011 ); Mugglinetal. ( 2000 ); Youngetal. ( 2009a ); Zhuetal. ( 2003 )),resultinginpoint-to-arealmisalignment.Misalignmentalsoarisesinothertypesofenvironmentalstudies.Forexample,whenrelatingsoilcharacteristicstorainfall,thesoilsamplesandrainfallmonitorsareobservedatdifferentpointsinspace( Zhangetal. ( 2011 )).Thedisparatedatasetsmustbealignedatacommonspatialsupportbeforeinferencecanbeconducted.Berksonerrorisinducedwhenusingkriging(orothersmoothingmethods)toalignsuchspatiallyreferenceddata( Gryparisetal. ( 2009 ); Lopianoetal. ( 2011 ); Madsenetal. ( 2008 ); Pacioreketal. ( 2009 ); Szpiroetal. ( 2011 )).Toillustrate,considerthefollowingargumentforthepoint-to-pointmisalignmentproblem(ananalogousargumentcaneasilybeshownforthepoint-to-arealmisalignmentscenariowhenblockkrigingisusedtopredictthecovariate'sarealunitaverage).SupposetheresponsevariableYisobservedatnspatiallocationssu,andtheexplanatoryvariableXisobservedatnodifferentspatiallocationsso.Supposethe 77

PAGE 78

relationshipbetweenY(su)andX(su)isthatgivenin 4 .LetX(su)betheunobservedXvaluesatthelocationssuandX(so)betheobservedXvaluesatthelocationsso.SupposethejointdistributionofX(su)andX(so)is 264X(su)X(so)375N0B@264C(su)C(so)375,X(x)=264X,uuX,uoX,ouX,oo3751CA,(4)whereC(su)andC(so)arenkandnokmatricesofobservedcovariatesatthelocationsinsuandso,respectively,isanunknownvectorofkparameters,xisanunknownvectoroflparametersthatdenethecovariancematrixXwhichispartitionedintothecovariancematricesX,uu,thecovarianceamongtheunobservedlocations,X,uo,thecovarianceamongtheunobservedandobservedlocations,andX,oo,thecovarianceamongtheobservedlocations.Forexample,xmaybeparametersthatdeneanexponentialcovariancestructure.Theconditional(kriging)distributionofX(su)jX(so)is X(su)jX(so),,xN(W=C(su)+X,uo)]TJ /F5 7.97 Tf 6.59 0 Td[(1X,oo(X(so))]TJ /F10 11.955 Tf 11.96 0 Td[(C(so))x=X,uu)]TJ /F10 11.955 Tf 11.96 0 Td[(X,uo)]TJ /F5 7.97 Tf 6.59 0 Td[(1X,ooX,ou). (4) NoticethekrigingmeanhasaBerksonerrorrelationshipwiththeunobservedcovariate, X(su)jX(so),,x=W+x,(4)wherexN(0,x).BecauseofthisBerksonerrorrelationship,Y(su)ExponentialFamily(mean=,scale=)g()=01n1+1W+1xxN(0,x). (4)Theframeworkhereisidenticaltotheframeworkdescribedin 4 withBerksonparameters=f,xg. 78

PAGE 79

Berksonerroringeneralizedlinearmodelshasbeenstudiedextensively.Regressioncalibration,simulation-extrapolation,quasi-likelihoodandvariancefunction(QVF),pseudo-maximumlikelihoodestimators(PMLEs)andmodiedestimatingequations(MEEs)( Buonaccorsi ( 2009 ); Carrolletal. ( 2006 ))havebeenconsideredforthisproblem.However,becausetheseapproachesgenerallyassumetheoff-diagonalelementsofxarezero,themethodsarenotimmediatelyapplicablefortheBerksonerrormodelthatarisesinspatialmisalignmentproblems.AlthoughtheBerksonerrorproblemhasbeenconsideredinthecontextofnormallinearmodelswithmisalignedspatialdata( Gryparisetal. ( 2009 ); Madsenetal. ( 2008 ); Szpiroetal. ( 2011 ),andChapter 3 ),Berksonerrorinducedfromaligningdisparatespatialdatasetsingeneralizedlinearmodelshasreceivedlessattention.Afullymodel-basedBayesianapproachwaspresentedin Mugglinetal. ( 2000 ). Gryparisetal. ( 2009 )brieydiscussedmeasurementerrorinprobitregressionmodelsandconductedasimulationstudyforalogisticregressionmodel.Inthesimulationstudy,theauthorscomparedthenaiveestimatorfoundbyusingthekrigingmeanasifitwerethetruevalue,anout-of-sampleregressioncalibration(RC-OOS)estimator( Thurstonetal. ( 2003 )),andanexposuresimulationapproach.Becausetheexposuresimulationapproachactuallyinducesbiasandresultsinanincorrectimputationscheme,theauthorsrecommendagainstthatapproach.Inthiswork,weconsiderthemodelin 4 wherethegoalistoconductinferenceabouttherelationshipbetweenYandXwhenonlyWisobserved.FornormallyandPoisson-distributedresponsesusingthecanonicallink,naiveanalysesresultinunbiasedestimatesoftheparameters.However,thestandarderrorsassociatedwiththeparameterestimateswillbeunderestimatediftheBerksonerrorisignored.Forageneralizedlinearmodelwithbinomiallydistributedresponsesusingthecanonicallink,naiveestimatesarebiasedandunderestimatetheuncertainty( Buonaccorsi ( 2009 )). 79

PAGE 80

WeillustratehowadditionaluncertaintyinducedfromusingWasthepredictorcanbeaccountedforingeneralizedlinearmodelsusingstandardtechniquesdevelopedformixedmodels.Indoingso,theparameterandanassociatedmeasureofuncertaintyareestimatedinasingleprocedure.Themethodsarevalidforhomoscedastic,heteroscedastic,andcorrelatedmeasurementerrors.Forconcreteness,weconsiderthespatialmisalignmentproblemasanexampleofaBerksonerrormodelinpracticeandillustratetheproposedmethodologywithtwoexamples:1)relatingFloridacountyasthmacountstoPM2.5concentrations,theconcentrationofparticlessmallerthan2.5micrometers;and2)relatingtheprobabilityoflowbirthweighttoPM2.5concentrations.Theperformanceoftheproposedmethodologyisassessedusingasimulationstudy.Thepaperisorganizedasfollows.InSection 4.2 ,theproposedmethodologyisdescribed.InSection 4.3 ,thenaivemethodthatignorestheBerksonerroriscomparedtotheproposedmethodologies.Then,asimulationstudyisdescribedinSection 4.4 .Finally,theapplicationsarepresentedinSection 4.5 followedbyadiscussioninSection 4.6 4.2ProposedMethodologyUndermodel 4 ,ourgoalistodevelopaframeworktoestimatetheparameter1andanappropriatemeasureofuncertaintythataccountsfortheadditionalerrortermx.Notethattherelationshipin 4 isidenticaltothespecicationofageneralizedlinearmixedmodelwherexplaystheroleofarandomeffect.Thusweproposeadaptingthemethodsusedtoestimateparametersingeneralizedlinearmixedmodelstoestimate.Whentheresponsesarenormallydistributedandanidentitylinkisused,themodelbecomes: YN(~,)~=01n1+1W=21x+y. (4) 80

PAGE 81

InChapter 3 ,weemployedmethod-of-momentsandrestrictedmaximumlikelihood-likemethodstoobtainanestimateofandthesubsequentlydenedestimatedgeneralizedsquaresestimator.Fornon-normalresponses,weconsideradaptingthepenalizedquasi-likelihood(PQL)methodof Breslow&Clayton ( 1993 ). Breslow&Clayton ( 1993 )proposedaFisherscoringalgorithmtoestimateregressionparametersandtopredictrandomeffectsinageneralizedlinearmixedmodel.ThealgorithmmustbeadaptedintwowaysforusewiththeBerksonerrormodelin 4 :1)Thepresenceof1inboththexedandrandomcomponentmustbeaccountedfor,and2)KnownorexternallyestimatedBerksonparametersmustbeincorporated.Thus,thealgorithmbecomes: 1. Initializeandxusingb(i)andb(i)x,thenaiveestimatesofandx,respectively 2. DenetheworkingresponseZ(i+1)=b(i)01n1+b(i)1W+b(i)1b(i)x+g0g)]TJ /F5 7.97 Tf 6.58 0 Td[(1b(i)01n1+b(i)1W+b(i)1b(i)xY)]TJ /F4 11.955 Tf 11.95 0 Td[(g)]TJ /F5 7.97 Tf 6.58 0 Td[(1b(i)01n1+b(i)1W+b(i)1b(i)xwhereg0isthederivativeofthelinkfunction.UnlikeintypicalGLMMswhereaknownmatrixismultipliedbytherandomeffect,1ismultipliedbytherandomeffect.Asaresult,hereweplug-intheworkingestimateof1. 3. EstimateD,theusualdiagonalGLMweightmatrix,asbD(i+1),wherebD(i+1)jj=wj Vg)]TJ /F5 7.97 Tf 6.58 0 Td[(1b(i)0+b(i)1Wj+b(i)1b(i)j,XjWg0g)]TJ /F5 7.97 Tf 6.58 0 Td[(1b(i)0+b(i)1Wj+b(i)1b(i)j,XjW2whereWjisthejthelementofthevectorW,wjarepotentialweightsandV()istheusualvariancefunction(see McCullagh&Nelder ( 1989 )). 4. Calculateb(i+1)=(bD(i+1)))]TJ /F5 7.97 Tf 6.59 0 Td[(1+b1(i)2x.IntypicalGLMMs,theparametersthatdenexareestimatedusingtheapproximateprolequasi-likelihoodfunction.Here,becausetheparametersthatdenexareassumedknownorestimatedexternally,theadditionalsteprequiredtoestimatetheparametersisomitted. 5. Calculateb(i+1)x=b(i)1x(b(i+1)))]TJ /F5 7.97 Tf 6.59 0 Td[(1(Z)]TJ /F16 11.955 Tf 11.95 0 Td[(Wb(i)),whereW=[1n1W] 6. Calculateb(i+1)=(W0(b(i+1)))]TJ /F5 7.97 Tf 6.59 0 Td[(1W))]TJ /F5 7.97 Tf 6.59 0 Td[(1W0(b(i+1)))]TJ /F5 7.97 Tf 6.58 0 Td[(1Z. 81

PAGE 82

7. RepeatSteps2-6untilconvergence.Thecovarianceofthenalestimatorisestimatedvia dVar(b)=(W0b)]TJ /F5 7.97 Tf 6.59 0 Td[(1W))]TJ /F5 7.97 Tf 6.58 0 Td[(1.(4)ForthemethodologiesinChapter 3 andhere,^21isusedtoapproximate21.Dependingontherelationshipamong^1,itsvariance,andx,anadjustmentmaybeneededtoaccountforthefactthatwhen^1isanunbiasedestimateof1,^21isabiasedestimateof21.Thus,ratherthanusing^21toestimate21,wecanuse ^21)]TJ /F14 11.955 Tf 12.1 3.16 Td[(dVar(^1)(4)because Var(^1)=E(^21))]TJ /F7 11.955 Tf 11.95 0 Td[(21.(4)However,if^21islargerelativetoitsvarianceestimate,theeffectwillbenegligible.Thealgorithmaboveassumesisknown,whichisrarely,ifever,thecase.However,ifadatasetthatdependsonisavailable,thenthedatasetcanbeusedtoobtainanestimateof.Thentheestimates^canbeusedasiftheyarethetruevaluesinthealgorithmdetailedabove.Historically,estimatorsfoundusinglikelihoodmethodswithestimatedparameterstreatedasknownarereferredtoaspseudo-likelihoodestimators( Gong&Samaniego ( 1981 )).Thus,werefertothethetwo-stepestimatordenedhereasapseudo-penalizedquasi-likelihoodestimator.WhenBerksonerrorarisesfromspatialmisalignment,theobservedX(so)canbeusedtoestimateandxwithEGLSandrestrictedmaximumlikelihood.Thatis,^xisfoundbymaximizingtherestrictedlikelihoodofX(so)and ^EGLS=(C(so)0^)]TJ /F5 7.97 Tf 6.58 0 Td[(1X,ooC(so)))]TJ /F5 7.97 Tf 6.59 0 Td[(1C(so)0^)]TJ /F5 7.97 Tf 6.59 0 Td[(1X,ooX(so),(4)where^X,oo=X,oo(^x).Then,thealgorithmabovecanbeusedwith^=f^EGLS,^xg.Here,asinChapter 3 ,weaccountfortheerrorassociatedwithestimating. 82

PAGE 83

Recall, 264X(su)X(so)375N0B@264C(su)C(so)375,X(x)=264X,uuX,uoX,ouX,oo3751CA,(4)and X(su)jX(so),,xN(W=C(su)+X,uo)]TJ /F5 7.97 Tf 6.59 0 Td[(1X,oo(X(so))]TJ /F10 11.955 Tf 11.96 0 Td[(C(so))x=X,uu)]TJ /F10 11.955 Tf 11.96 0 Td[(X,uo)]TJ /F5 7.97 Tf 6.59 0 Td[(1X,ooX,ou). (4) Toremoveconditioningon,considerauniformpriorfor.TheresultingposteriordistributionjX(so),xis jX(so),xN(~GLS=(C(so)0)]TJ /F5 7.97 Tf 6.58 0 Td[(1X,ooC(so)))]TJ /F5 7.97 Tf 6.59 0 Td[(1C(so)0)]TJ /F5 7.97 Tf 6.59 0 Td[(1X,ooX(so),(C(so)0)]TJ /F5 7.97 Tf 6.58 0 Td[(1X,ooC(so)))]TJ /F5 7.97 Tf 6.59 0 Td[(1). (4) Theconditionaldensityis (X(su)jX(so),x)=Z(X(su)jX(so),x,)(jX(so),x)d, (4) wherethenotation(xjy)representsthedensityofxgiveny.Basedontheintegral 4 X(su)jX(so),xN(~W=C(su)~GLS+X,uo)]TJ /F5 7.97 Tf 6.59 0 Td[(1X,oo(X(so))]TJ /F10 11.955 Tf 11.95 0 Td[(C(so)~GLS),~x=X,uu)]TJ /F10 11.955 Tf 11.96 0 Td[(X,uo)]TJ /F5 7.97 Tf 6.59 0 Td[(1X,ooX,ou+(C(su))]TJ /F10 11.955 Tf 11.95 0 Td[(X,uo)]TJ /F5 7.97 Tf 6.59 0 Td[(1X,ooC(so))(C(so)0ooC(so)))]TJ /F5 7.97 Tf 6.58 0 Td[(1(C(su))]TJ /F10 11.955 Tf 11.95 0 Td[(X,uo)]TJ /F5 7.97 Tf 6.59 0 Td[(1X,ooC(so))0). (4) NowtheBerksonparametersarexandthealgorithmcanbereformulatedbyreplacingWandxwith~Wand~x,respectively.Indoingso,theestimatorsaccountfortheclassicalmeasurementerrorassociatedwithestimatingbutnotthatassociatedwithestimatingx.Althoughitwouldseemreasonabletoremovethedependence 83

PAGE 84

onxaswell,theposteriordensitydoesnotexistinclosedformandafullyBayesianapproachwouldbemoreappropriate.ForadiscussionofBayesianmethodsinspatialmisalignmentproblemssee Gelfandetal. ( 2001 ); Mugglinetal. ( 2000 ); Zhuetal. ( 2003 ). 4.3ComparingNaiveandAdjustedEstimatorsHereweconsidertheeffectofusinganaiveestimator.Naiveestimatorsaretheestimatorscalculatedwhenoneignoresthemeasurementerror,usesthesurrogateWasifitwerethetruevalue,andconductsinferenceassumingtheoriginalregressionmodel.Fornormallydistributedresponses,thenaiveestimatoroftheregressionparametersisthegeneralizedleastsquaresestimator ^naive=(W0)]TJ /F5 7.97 Tf 6.58 0 Td[(1yW))]TJ /F5 7.97 Tf 6.59 0 Td[(1W0)]TJ /F5 7.97 Tf 6.58 0 Td[(1yY.(4)Thevariance^naiveistakentobe(W0)]TJ /F5 7.97 Tf 6.59 0 Td[(1yW))]TJ /F5 7.97 Tf 6.59 0 Td[(1.However,becausethevarianceestimatorisbasedontheoriginalvarianceofY,y,thisisnotcorrect.Thevarianceof 4 is Var(^naive)=(W0)]TJ /F5 7.97 Tf 6.59 0 Td[(1yW))]TJ /F5 7.97 Tf 6.59 0 Td[(1W0)]TJ /F5 7.97 Tf 6.59 0 Td[(1y)]TJ /F5 7.97 Tf 6.59 0 Td[(1yW(W0)]TJ /F5 7.97 Tf 6.58 0 Td[(1yW))]TJ /F5 7.97 Tf 6.59 0 Td[(1.(4)Thisisapproximatelyequalto (W0)]TJ /F5 7.97 Tf 6.58 0 Td[(1yW))]TJ /F5 7.97 Tf 6.59 0 Td[(1(4)if )]TJ /F5 7.97 Tf 6.59 0 Td[(1y=)]TJ /F5 7.97 Tf 6.59 0 Td[(1y(21x+y)I,(4)orequivalently,if(21x+y)y.Bothyandmustbeknowntoassess 4 .However,ifbothmatricesareknown,thenthecorrectGLSestimatorcouldbecalculated.Ifthestructuresofyandaresimilar,thenthenaiveapproachwouldbeapproximatelycorrectasdescribedinmoredetailnext. 84

PAGE 85

CaseI:y=2yIandx=2xITheOLSestimatorwouldbethepointestimatorforboththenaiveandadjustedmethod.However,thenaivevarianceestimator, 2y(W0W))]TJ /F5 7.97 Tf 6.59 0 Td[(1,(4)wouldunderestimatethetruevariability, (212x+2y)(W0W))]TJ /F5 7.97 Tf 6.59 0 Td[(1.(4)Theratio r=2y 212x+2y(4)quantiestheproportionbywhichthevariabilityisunderestimated.Ifr1then212x<<2yand2y212x+2y.Therefore,2yIdominatestheinducederrorstructure.Inthiscase,using 4 wouldnotbeproblematicwhenconductinginference.Further,because2xand2yareunknownandboth y=2yI(4)and =(212x+2y)I(4)havethesamestructure,aconstanttimesadiagonalmatrix,themeansquarederror(MSE)willbeunbiasedfor(212x+2y).Thus,regardlessofr,theestimatedvariancewillbeunbiasedforthetruevarianceoftheestimatorandsubsequentinferenceswillbevalid.Thatis,becausethestructuresofthetwocovariancematricesarethesame,theinferenceproceduresperformasdesired,andtheBerksonerrorcanbeignoredinpractice.CaseII:GeneralyandxForthecaseofgeneralcovariancematricesyandx,ifthestructureofyandaresimilar,thennaiveapproachwouldbeapproximatelycorrect.Asanexample, 85

PAGE 86

supposey=2yR(y)andx=2R(x).Ifx=y=,thentheinducedcovariancematrixbecomes =(212x+2y)R().(4)Becauseandyhavethesamestructure,aconstanttimesthematrixR(),naiveinferencebasedonassumingthisparametricstructureforthecovariancematrixwillbevalid.However,thesimilarityofthetwomatriceswillbeunknowninpractice.Thus,itisimpossibletodeterminepriortoconductinginferencewhethertheBerksonerrorcanbeignored.CaseIII:SpatialMisalignmentwithindependentregressionerrorsConsideraspatialmisalignmentscenariowithy=2Iandx=X,uu)]TJ /F10 11.955 Tf 11.95 0 Td[(X,uo)]TJ /F5 7.97 Tf 6.58 0 Td[(1X,ooX,ou.Thecovariancematrixis=21x+2I.NaivemethodswillbevalidifcIwherecisaconstant.Becauseisapositive-denitecovariancematrix,alleigenvaluesofarerealandpositive.Inaddition,alleigenvaluesofarecontainedinatleastoneGershgorindisc,[i,i)]TJ /F4 11.955 Tf 12.5 0 Td[(Ri,i,i+Ri],i=1,...,n,whereRi=Pj6=iji,jjandi,jisthe(i,j)thelementofthematrix(seeTheorem7.2.1 Golub&vanLoan ( 1996 ),p.320).Infact,itcanbeshowntheeigenvaluesofareboundedbetween[2,2+Rmax]whereRmax=maxRi(seeAppendix).Whenalleigenvaluesareequal,thecovariancematrixisaconstanttimestheidentity.Therefore,asRmaxincreases,theeigenvaluesbecomemoredisparateandthematrixbecomeslessdiagonalinthesenseitcannotbeapproximatedwellbyaconstanttimestheidentity.BecauseRi=Pj6=iji,jj,therangeoftheobservedeigenvalueswillincreaseastheoffdiagonalelementsofincreaseinmagnitudeandasthenumberofnon-zerooff-diagonalelementsincreases.IfthecovariancestructureamongtheX 86

PAGE 87

isanexponentialonewithscale2xandrangex,then i,j=212x[exp()]TJ /F4 11.955 Tf 9.3 0 Td[(di,j=x))]TJ /F11 7.97 Tf 16.28 14.95 Td[(noXk=1noXl=1exp()]TJ /F6 11.955 Tf 9.3 0 Td[((di,l+dj,k)=x)D)]TJ /F5 7.97 Tf 6.59 0 Td[(1oo(i,j)],(4)whereDoo(l,k)=exp()]TJ /F4 11.955 Tf 9.3 0 Td[(dl,k=x)anddi,jistheEuclideandistancebetweenlocationsiandj.Thus,as1and2xincreaseinmagnitude,theabsolutevalueofi,jwillincrease.Forxedno,ifnu,thenumberoflocationswheretheresponseisobserved,increases,thenRiwillincreaseforalli.Themagnitudeoftheincreasedependsonthenumberoflocationswherethecovariateisobservedandthedistancebetweenthelocationswheretheresponseisobservedandthelocationswherethecovariateisobserved.IncreasingbothnuandnomayresultineitheranincreaseordecreaseintheRidependingonthelocationswheretheobservationsareadded.Asaconsequence,increasingthenumberofobservationswhereYisobserved,increasing2x,orincreasingj1jwillincreasetheeffectofBerksonerrorandcauseestimatesbasedonthenaivemethodstoperformpoorly(e.g.,biasedestimatesofthestandardandcoveragebelownominallevels).PoissonandBinomialFornon-normallydistributedresponses,thecomparisonofnaiveinferenceandPQL-basedinferencerequiresadifferentapproach.IgnoringtheBerksonerror,thecovarianceoftheestimatedregressionparametersis Var(^naive)=(W0DW))]TJ /F5 7.97 Tf 6.59 0 Td[(1,(4)whereDisadiagonalmatrixoftheGLMweights.InthecaseofPoissonregressionwithaloglink,theithdiagonalelementofDisexp(0+1Xi),andinthecaseoflogisticregression,theithdiagonalelementofDisnTexp(0+1Xi) (1+exp(0+1Xi))2, 87

PAGE 88

wherenTisthenumberoftrials.ThePQL-basedvarianceestimateis Var(^PQL)=(W0)]TJ /F5 7.97 Tf 6.58 0 Td[(1W))]TJ /F5 7.97 Tf 6.59 0 Td[(1,(4)where=D)]TJ /F5 7.97 Tf 6.59 0 Td[(1+21x.Theratioofthevariancesof^1is r=(W0DW))]TJ /F5 7.97 Tf 6.59 0 Td[(1[2,2] (W0)]TJ /F5 7.97 Tf 6.58 0 Td[(1W))]TJ /F5 7.97 Tf 6.59 0 Td[(1[2,2],(4)wherethesubscript[i,j]referstothe(i,j)thelementofthematrix.Because 4 isdifculttocomputeinclosedform,weconsideranapproximationasfollows.SupposeDc1Iwherec1isaconstant.Further,supposec2Iwherec2isaconstant.ForthecaseofPoissonregression,becausetheXivary,weusetheexpectedvalueofXitodenec1=exp(0+1E(Xi)).Inthecontextofspatialmisalignment,ifthecovariancematrixforthecovariateisdenedbyanexponentialcovariancematrixwithscale2xandrangex,thenwedenec2=1 c1+212x1)]TJ /F6 11.955 Tf 11.96 0 Td[(exp)]TJ /F4 11.955 Tf 13.09 8.08 Td[(k x,wherekisanunknownvalue.Theterm2x1)]TJ /F6 11.955 Tf 11.96 0 Td[(exp)]TJ /F4 11.955 Tf 13.09 8.09 Td[(k xisusedtoapproximateadiagonalelementofx. 88

PAGE 89

Basedontheapproximations, rc)]TJ /F5 7.97 Tf 6.59 0 Td[(11 1 c1+212x(1)]TJ /F6 11.955 Tf 11.96 0 Td[(exp()]TJ /F11 7.97 Tf 12.6 4.7 Td[(k x)) (4) =c11 c1+212x1)]TJ /F6 11.955 Tf 11.95 0 Td[(exp)]TJ /F4 11.955 Tf 13.09 8.08 Td[(k x)]TJ /F5 7.97 Tf 6.59 0 Td[(1 (4) =1 1+c1212x1)]TJ /F6 11.955 Tf 11.95 0 Td[(exp)]TJ /F11 7.97 Tf 12.59 4.7 Td[(k x (4) approximatestheratioofthenaivevarianceestimatetothePQLvarianceestimate.Forthelogisticregressionmodel,thesamerelationshipholdswithc1=nTexp(0+1E(Xi)) (1+exp(0+1E(Xi)))2.Inthenextsection,asimulationstudyisconductedtoincreaseourunderstandingofhownaiveandadjustedestimatorsperforminpracticeandassesstheapproximationsdescribedhere. 4.4SimulationStudyWeconsiderapoint-to-pointmisalignmentprobleminthesimulationstudy.Thespatialdomainis[0,500][0,500]R2.Theresponseisobservedat100locationsandthecovariateisobservedat100differentlocations.ThesurfaceXwasgeneratedfromX(s)GaussianProcess(=40,X),whereXisdenedbyanexponentialcovariancestructureparameterizedbyascaleandrangeparameter.Weconsiderthreedifferentdistributionsfortheresponseandusethecanonicallink:1)Poisson(loglink),2)Binomial(logitlink),and3)Normal(identitylink).ForPoissonandlogisticregression,theusualGLMweightsdependonthemeanparameters0and1.TheweightsderivedusingthePQLalgorithmdependon0,1,thescaleandtherange.Forthelinearmodel,weconsideredaniidmodelfortheoriginalregressionerror,N(0,2I).Foralldistributionswevary1,thescale,therangeandtheBerksonparameters: 1:4levels:0.02,0.05,0.07,0.10; 89

PAGE 90

Scale:PoissonandBinomial-2levels:35,70Normal-3levels:35,70,105; Range:3levels:50,200,400; BerksonParameters:2levels:knownandunknown.Forthenormaldistribution,theparameter0wasxedto0and4valuesof2wereconsidered,2,4,8,and16.ForthePoissondistribution,14levelsof0wereconsideredandfortheBinomialdistribution,8levelsof0wereconsidered: Poisson:denedas0=log(M))]TJ /F7 11.955 Tf 12.28 0 Td[(140:14levelsofM:1,2,4,8,12,16,20,24,28,32,36,40,44,48; Binomial:denedas0=log(M) 1)]TJ /F5 7.97 Tf 6.59 0 Td[(log(M))]TJ /F7 11.955 Tf 12.21 0 Td[(140:8levelsofM:0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8.Inaddition,forthebinomiallydistributedresponseweconsideredthreelevelsofthenumberoftrialsateachlocation:1,10,50.Crossingthefactorsleadstoatotalof423214=672settingsforthePoissondistribution,423283=1152settingsfortheBinomialdistribution,and43324=288settingsforthenormaldistribution.Foreachofthe2112simulationsettings,1000simulateddatasetsweregenerated.Asdiscussedlater,additionalnormalsimulationswerealsoconducted. 4.4.1StudyResultsBecausetheresultsarequalitativelysimilarforknownandunknownBerksonparameters,allguresandtablescorrespondtothesettingswheretheBerksonparameterswereknown.Asubsetofresultsarepresentedhere,buttheremainingsimulationresultssupporttheconclusionsdrawn.Ineachtable,themean(^1),variance(^2^1),bias,andmeansquarederrorassociatedwiththe1000estimatesof1areprovided.Inaddition,theaverageoftheestimatedvariancesof^1(s2^1)andtheproportionofthe95%condenceintervalscoveringthetruevalueareprovided. 90

PAGE 91

4.4.1.1PoissonForboththenaivemethodandthePQLmethod,theestimatesoftheregressionparameter1werevirtuallyindistinguishable,regardlessofthevaluesoftheotherfactors.Forbothmethods,thebiaseswerenegligible.Forthenaivemethod,thevarianceestimatorunderestimatedthetruevariability.Theeffectofeachfactorgivenallotherfactorsxedwasexplored.ThemagnitudeofthebiasincreasedastheexpectedvalueofYatX=40increased,astheslopeincreased,astherangedecreased,orasthescaleincreased(seeTable 4-1 ).Inaddition,thevariabilityoftheestimatorwaslargerwhentheBerksonparameterswereestimated.Becausethenaiveestimateofthevariabilitywasdownwardlybiased,coveragefor95%condenceintervalswerebelow95%andfollowedthesamedecreasingrelationshipasthebiasinthevarianceestimator(Figure 4-1 ).ForthePQLmethod,thevarianceestimatorapproximatedthetruevariabilitywell.ThevariabilitywaslargerwhentheBerksonparameterswereestimated.Biaseswerenegligibleandcoveragewasnear95%inallscenarios(Figure 4-2 andTable 4-1 ).WhencomparingthenaivevarianceestimatestothePQLvarianceestimates,bothwhentheBerksonparameterswereknownandunknown,theratiooftheestimatesofvarianceof^1followedtherelationshipderivedinSection 4.3 .Here,wechosek=40byvisuallyassessingthet.AsshowninFigure 4-3 ,thecurvesapproximatedtherelationshipwell,butnotperfectly. 4.4.1.2BinomialAlthoughnotdisplayedhere,for1trial,thenaiveandPQLestimatesoftheslopeanditsvariancewerenotdifferent.Biasesinboththeparameterestimateanditsvariancedidnotaffecttheperformanceof95%condenceintervals.Thatis,coveragefor95%condenceintervalswasnear95%forallsettings(resultsnotshown).For10trials,inalmostallsettings,thenaiveandPQLestimatesoftheslopeanditsvariancewerenotdifferentandbiasesdidnotaffectcoverage(seeTable 4-2 ).Theonly 91

PAGE 92

exceptionwasthecasewhentheslopewasequalto0.1andthescalewasequalto70.Inthiscase,thenaivevarianceestimatewasdownwardlybiasedandcoveragefor95%intervalswerebetween80%and90%.ThePQLmethodperformedwellinthiscase,andcoverageswerenear95%.With50trials,forthenaivemethod,thebiasinestimatesof1werenegligible.However,theassociatedestimatesofthevarianceunderestimatedthetruevariability.WhentheexpectedvalueofYatX=40wasequalto0.5,thebiasinthevariabilityreachedamaximumanddecreasedastheexpectedvalueofYatX=40movedawayfrom0.5.Thiswasexpectedbasedonthesymmetryofthevariancefunctionaround0.5.Withrespecttotheothervariables,thesamerelationshipsdiscussedforthenaivevarianceinthePoissonsettingswereobserved.Asaresult,coveragefor95%condenceintervalsisbelow95%(seeTable 4-2 ).ThePQL-basedestimatesof1wereessentiallyunbiased.Estimatesofthevariabilityweredownwardlybiasedwhentherangewas50andupwardlybiasedwhentherangewas200or400.ThebiasdecreasedastheexpectedvalueofYatX=40increased.However,becausethebiaswasrelativelysmall,coveragefor95%condenceintervalswasnear95%forallsettings.IncomparingthenaivevarianceestimatestothePQLvarianceestimates,bothwhentheBerksonparameterswereknownandunknown,theratiooftheestimatesofvarianceof^1followedtherelationshipderivedinSection 4.3 .Again,wechosek=40basedonvisuallyinspectingthet.AsshowninFigure 4-4 ,thecurvesapproximatedtherelationshipwell,butnotperfectly. 4.4.1.3NormalFortheparametervaluesconsidered,themaximumradiioftheGershgorindiskswaslessthan2forallpossiblesettings.Asaresult,theinducedcovariancematrixwasapproximatelydiagonalandnaivemethodswereexpectedtobeperformwell. 92

PAGE 93

Thesimulationresults(notshown)conrmednaivemethodsperformedwell,andweconsideredadditionalsettingstoillustratewhennaivemethodsareinadequate.Wexed0=0,2=16,thescaleequalto105,andtherangeequalto200.WeincreasedthenumberoflocationswhereYwasobservedto300.Toillustratetheeffectof1onnaivecoverage,weconsidered1equalto0.1,0.2,0.4,0.6and0.8.ThenaivecoveragewithknownBerksonparametersandthecorrespondingmaximumGershgorinradiusareshowninTable 4-3 .Asexpected,as1increases,themaximumGershgorinradiusincreasesandnaiveinferenceperformsmorepoorly.Weprovidemoredetailforthecase1=0.8inTable 4-4 .Thenaiveestimatesofthevarianceunderestimatedthetruevariability.Asaresult,coveragefor95%condenceintervalswaslessthan95%.FortheMOMandREML-likemethods,estimatesofthevariabilitywereessentiallyunbiasedandcoveragefor95%condenceintervalswerenear95%inallcases. 4.5ApplicationIn2007,duringthemonthsofApril,MayandJune,thelargestreinthehistoryofFloridaandGeorgia,theBugabooScrubFire,ragedinsouthernGeorgiaandnorthcentralFlorida.SmokefromtherereachedAtlanta,southFloridaandevenpartsofAlabamaandMississippi( GeorgiaForestryCommision ( 2007 )).SmokefromtheearlystagesoftherepouredintoFloridafromGeorgia.OnApril29th,thesmokewasrelativelyconcentrated.However,onApril30th,smokefromthereblanketedmostofnorthcentralFloridaandevenpartsofcentralFlorida.OnMay3rd,thenGovernorCharlieCristsignedexecutiveorder07-86promptingareassistancegrantfromtheFederalEmergencyManagementAgency( TheStateofFloridaStateEmergencyResponseTeam ( 2007 )).OnMay10th,astatewidehealthadvisorywasreleasedbytheFloridaDepartmentofEnvironmentalProtection:CautionStatement:Theelderly,children,andthosewithrespiratorydiseasesincludingasthmaorheart 93

PAGE 94

diseaseshouldavoidprolongedoutdooractivity;everyoneelseshouldlimitprolongedoutdooractivity.( FloridaDepartmentofEnvironmentalProtection ( 2007 )) 4.5.1Point-to-ArealMisalignmentInlightofthesuddenonsetofdrasticairqualityconditionsattheendofAprilandtheincreasedpublicawarenessofthewildresbyMay10th,wehypothesizedasthmahospitalandemergencyroomadmissionswouldbecorrelatedwithPM2.5duringtheearlystagesofthere.Becausebehaviorlikelychangedshortlyaftertheonsetoftherstextremesmokeeventandtheemergencydeclaration,weconsideredthetotalnumberofhospitalizationsandERcasesbetweenApril30th,2007andMay6th,2007.ThePMvaluedenedasthemaximumvalueofPM2.5between4/30/2007-5/6/2007wasobservedateachof32airqualitymonitorsinFlorida(Figure 4-5 )( U.S.EnvironmentalProtectionAgency ( 2007 )).Florida'sDepartmentofHealth(FDOH)hasadata-sharingagreement(DSA)withFlorida'sAgencyforHealthCareAdministration(AHCA),whichprovidedaccesstotwodatasources:condentialhospitalizationrecordsandemergencyroomrecords.Thehealthdataconsistedofallrecordswhereeithertheprincipalorotherdiagnosiscodeswereforasthma,ICD-9codes493.0-493.9or786.07( NationalCenterforHealthStatistics ( 1998 )).Toaccountfordemographicvariation,theexpectednumberofhealtheventswerecalculatedforeachcountybycalculatingthestandardizedeventratioadjustingforage,genderandethnicityusingthestateofFloridaasthestandardpopulation(seeAppendix).LetY(c)andE(c)denotetheobservedandexpectednumberofcasesincountyc,respectively,c=1,...,67.LetX(s)denotethePMvalueatlocations,s2DR2,whereDdenotesthestateofFlorida.Suppose X(s)GaussianProcess(,X(x)),(4) 94

PAGE 95

wherexisavectoroftheparametersthatdeneanexponentialcovariancestructure,i.e.,thescaleandrange.Foreachcountyc X(c)=jcj)]TJ /F5 7.97 Tf 6.59 0 Td[(1Zs2cX(s)ds(4)istheaveragePMvalueincountyc,wherejcjistheareaofcountyi.Letsodenotethecoordinates(latitude/longitude)ofthen=32monitors.LetthevectorX(so)bethevectorofthePMvaluesatthelocationsso.Note 4 implies X(so)N(1321,X,oo(x)).(4)Theparametersandxwereestimatedbasedonrestrictedmaximumlikelihoodandestimatedgeneralizedleastsquaresas^=27.89and^x=[186.72,1.42].LetCcorrespondtothe67countiesinFlorida.LetX(C)=[X(1),X(2),...,X(67)]0.Theblockkrigingdistributionis X(C)jX(so),,xN(X(C)=1671+C,so)]TJ /F5 7.97 Tf 6.59 0 Td[(1X,oo(X(so))]TJ /F10 11.955 Tf 11.96 0 Td[(1n),x=C,C)]TJ /F10 11.955 Tf 11.96 0 Td[(C,so)]TJ /F5 7.97 Tf 6.59 0 Td[(1X,ooso,C), (4) whereC,Cisthecovarianceamongtheblocks,C,soisthepoint-to-blockcovarianceand,X,ooisthecovariancematrixamongthepoints.Thedistributionin 4 canbewrittenas X(C)jX(so),,x=X(C)+v,vN(0,x).(4)ToapproximateX(C),rst,anegridwaslaidoverstate.Letsudenotethecoordinatesofthenu=3061gridlocations.Then,theconditionalmean, E(X(su)jX(so))=1nu+X,uo)]TJ /F5 7.97 Tf 6.58 0 Td[(1X,oo(X(so))]TJ /F10 11.955 Tf 11.95 0 Td[(1n),(4)wasestimatedusingtheestimatedvaluesofandx.Finally,X(C)wasestimatedas^X(C)byaveragingtheestimatedconditionalmeanvalueswithineachcounty.The 95

PAGE 96

blockkrigingvariancematrixwascalculatedanalogously.See Schabenberger&Gotway ( 2005 )formoredetailsregardingblockkriging.WemodeledthenumberofcasesusingPoissonregressionwiththeexpectednumberofcasesasanoffsetandaneffectduetoPM.BecauseX(c)isnotobserved,weinsteadconsideredthemodel Y(C)Poisson() (4) log()=log(E(C))+01671+1X(C)+1v (4) vN(0,x). (4) Toestimatetheparametersinthemodel,weusedtheestimatesofX(C)andxasiftheywereknown,andusedthepenalizedquasi-likelihoodalgorithmfromSection 4.2 toaccountforadditionalterm1v.Inaddition,wecarriedouttheregressionasif^X(C)werethetruecovariateX(C).Werefertothisasthenaiveapproach.Forthenaiveapproach,theestimate1was0.019withavarianceandstandarderrorof5.3e-6and0.0023,respectively.Whenusingpenalizedquasi-likelihood,theestimateof1was0.018withavarianceandstandarderrorof8.5e-6and0.0029,respectively.Accountingfortheadditionaluncertaintyduetokrigingresultedina60%increaseinthevarianceoftheestimateof1,5.3e-6to8.5e-6,anda27%increaseinthestandarderror,0.0023to0.0029.Theparameterestimatedecreased8%from0.019to0.018.AWaldtestforthesignicanceoftheparameterresultedinaz-scoreof8.41inthenaivecaseand6.09inthePQLcase,a27%decrease.However,inbothcases,theeffectofPMwassignicant.BasedonthePQLestimateof1,duringthetimeperiodunderconsideration,aunitincreaseinPMresultsinanestimatedmultiplicativeeffectof1.018ontherateofasthmacaseswitha95%condenceintervalof(1.012,1.024). 96

PAGE 97

4.5.2Point-to-PointMisalignmentThetimeperiodconsistingofweeksthreethrougheightofhumangestationisknownasembryogenesis.Duringembryogenesis,vitalorgansandlimbstructuresdevelop.Attheendofembryogenesis,theoffspringisofciallyreferredtoasafetus.Pregnancyisconsideredfulltermafter37weeksofgestation.Duringembryogenesis,theembryoisvulnerabletodevelopmentaldefects.Becauseofthis,wehopedtoassesswhetherhighlevelsofparticulatematterduringembryogenesisarerelatedtoanincreasedprobabilityoflowbirthweightwherealowbirthweightisdenedasafetusbeingbornweighinglessthan2500grams.ThePMvaluedenedasthemaximumamountofPM2.5duringthemonthofMaywasobservedateachof32airqualitymonitorsinthestateofFlorida(Figure 4-5 )( U.S.EnvironmentalProtectionAgency ( 2007 )).In2007,infantsofgestationalage37weeksbornbetweenDecember24thandDecember31stexperiencedthemajorityofembryogenesisduringthemonthofMay.ThebirthweightdataobtainedfromtheFloridaDepartmentofHealth'sOfceofVitalStatisticsconsistedof179non-Hispanicwhitemotherswhogavebirthduringthistimeperiodtooneinfant37weeksingestationalage.LetX(s)denotethePMvalueatlocations,s2DR2,whereDdenotesthestateofFlorida.Suppose X(s)GaussianProcess(,X(x)).(4)wherexisavectoroftheparametersthatdeneanexponentialcovariancestructure,i.e.,thescaleandrange.Letsodenotethecoordinates(latitude/longitude)ofn=32themonitors.LetthevectorX(so)bethevectorofthePMvaluesatthelocationsso.Theparametersandxwereestimatedbasedonrestrictedmaximumlikelihoodandestimatedgeneralizedleastsquaresas^=53.28and^x=[856.04,0.36]. 97

PAGE 98

Letsudenotethecoordinatesofthenu=179motherresidentiallocations.Theconditional(kriging)distributionofX(su)jX(so)is X(su)jX(so),,xN(X(su)=11791+X,uo)]TJ /F5 7.97 Tf 6.59 0 Td[(1X,oo(X(so))]TJ /F10 11.955 Tf 11.96 0 Td[(1321),x=X,uu)]TJ /F10 11.955 Tf 11.96 0 Td[(X,uo)]TJ /F5 7.97 Tf 6.59 0 Td[(1X,ooX,ou). (4) LetY(si)beanindicatoroflow-birthweightfortheithinfant,i=1,...,179.SupposetherelationshipbetweenY(su)andX(su)followsalogisticregressionmodel.Othercovariates(mother'sage,mother'seducation,mother'ssmokingstatus,mother'salcoholuse,mother'sfoodstampstatus)wereconsidered,butwerenotsignicantata10%signicancelevel.ThemodeldescribedherehasonlyPMasacovariate: Y(su)Bernoulli((su)) (4) logit((su))=011791+1X(su) (4) BecauseX(su)isunobserved,weinsteadconsideredthemodel Y(su)Bernoulli((su)) (4) logit((su))=011791+1X(su)+1v, (4) vN(0,x). (4) Toestimatetheparametersinthemodel,weusedtheestimatesofX(su)andxasiftheywereknown,andusedthepenalizedquasi-likelihoodalgorithmfromSection 4.2 toaccountforadditionalterm1v.Inaddition,weconsideredthenaivemodelthattreatsX(su)asthetruecovariate.Forthenaiveapproach,theestimateof1was0.056withavarianceandstandarderrorof0.000308and0.018,respectively.Whenusingpenalizedquasi-likelihood,theestimateof1was0.058withavarianceandstandarderrorof0.000370and0.019,respectively. 98

PAGE 99

Accountingfortheadditionaluncertaintyduetokrigingresultedina20%increaseinthevarianceoftheestimateof1,0.000308to0.000370,anda10%increaseinthestandarderror,0.018to0.019.Theparameterestimateincreased3%from0.056to0.058.AWaldtestforthesignicanceoftheparameterresultedinaz-scoreof3.19inthenaivecaseand3.00inthePQLcase,a6%decrease.However,inbothcases,theeffectofPMwassignicant.UsingthePQLbasedestimateof1,duringthetimeperiodunderconsideration,weestimatedaunitincreaseinPMresultedinanestimatedmultiplicativeeffectontheoddsoflowbirthweightof1.059witha95%condenceintervalof(1.02,1.10). 4.6ConcludingRemarksForlinearmodelsandgeneralizedlinearmodels,ifacovariatehasaBerksonerrorrelationshipwithitssurrogate,thengeneralizedleastsquaresandpenalizedquasi-likelihoodcanbeaugmentedandusedtoconductinference.Onecaveat,however,isinpracticetheparametersthatdenethesurrogateorthecovarianceassociatedwiththeerrormodelaretypicallyunknown.AlimitationofusingestimatedBerksonparametersisestimatingintroducesacomponentofclassicalmeasurementerrorthatcancausebiasinparameterestimatesandstandarderrorstobetoosmall.IftheBerksonparameterscanbeestimatedwithsufcientprecision,thentheapproachpresentedhereshouldbeadequate.Ifnot,bootstrapalgorithmscanbeusedinconjunctionwithmethodpresentedinSection 4.2 inanattempttoaccountfortheadditionalerror.However,forcomplicateddata-generatingmechanisms,derivingabootstrapprocedureisnottrivial.Inthecontextofmeasurementerrormodels,see Buonaccorsi ( 2009 )Section16.6foradiscussionofbootstrapmethods.Inourwork,weonlyconsideredPQLtoaccountfortheadditionaluncertaintyinthegeneralizedlinearmodel.ThePQLapproachcanperformpoorly,particularlywhenestimatingvariancecomponents( Agresti ( 2002 ),Section12.6.4).However,becausethevariancecomponentsareestimatedusinganexternaldatasetandused 99

PAGE 100

inthepseudo-PQLFisherscoringalgorithmasiftheyareknown,thepseudo-PQLapproximationforestimatingthemaineffectsproposedheredoesnotusevariancecomponentsthatareestimatedviarestrictedpenalizedquasi-likelihoodmethods.Asaresult,asillustratedinoursimulationstudies,thepseudo-PQLmethodperformswellwhenconductinginference.Nonetheless,exploringtheperformanceofothermethodsdevelopedforGLMMswouldbeofinterest.Inthecontextofspatialmisalignment,becausetheunobservedcovariatehasaBerksonerrorrelationshipwiththekrigingmean,themethodscanbeappliedtomisaligneddatainlinearandgeneralizedlinearmodels.However,additionalerrorisinducedfromestimatingthekrigingmeanandkrigingvariance,i.e.,theBerksonparameters.Althoughweaccountedforestimatingthemeanparameters,theadditionaluncertaintyinestimatingxwasnotaccountedfor.Toaccountforthisinlinearmodels, Szpiroetal. ( 2011 )developedthreebootstrapmethodologies.Althoughbootstrapmethodscouldbedevelopedfornon-normalresponses,inthesimulationstudyconsideredhere,fortheproposedmethodologies,theadditionaluncertaintyfromestimatingtheBerksonparametersprovednegligible.ThesimulationresultsalsoillustratenaiveinferenceinPoissonregressionmodelswithBerksonerrorcanperformpoorlydependingonthesizeoftheeffect,thesizeoftheintercept,andthenatureoftheBerksonerror.Consequently,incorrectinferencesmaybedrawnwhenusingPoissonregressioninthepresenceofBerksonerror.TheresulthasimplicationsforenvironmentalassociationstudieswhererelationshipsamongmisaligneddataareroutinelyassessedusingmethodsthatinduceBerksonerror.Forbinomiallydistributedresponses,thesimulationresultsillustratethat,asthenumberoftrialsincreases,theeffectofBerksonerrorbecomesmorepronounced.Finally,fortheindependentnormallydistributedresponses,thesimulationresultsshowthatincreasingtheabsolutevalueoftheoffdiagonalelementsoftheinducedcovariancematrixresultsindecreasingprobabilitiesofcoveragefornominal95%condenceintervals.Forthe 100

PAGE 101

spatialmisalignmentproblem,theperformanceofnaivemethodsdependsontheparameters1,therange,thescale,thenumberoflocationswhereYandXisobservedandthedistancebetweentheselocations.Here,weillustratedcoveragefor95%condenceintervalsfellbelowthenominallevelwhenweincreasedthemagnitudeof1andthenumberoflocationswhereYisobserved.Inpractice,however,withoutknowingthetrueparameters,thetruecoverageofnaivecondenceintervalswillbeunknownandwerecommendusinganadjustedproceduretoconductinference.Finally,BayesianhierarchicalmodelsareanappealingapproachtoaccountforBerksonerrorbecausetheparameterestimationprocesswouldnaturallypropagatealltypesofuncertaintythroughthemodel.However,suchanapproachcanbecomputationallyintensiveand,incomplexmodels,substantialeffortmayberequiredtosetupalgorithmstoapproximateposteriordistributions.AlthoughMarkovChainMonteCarlo(MCMC)methodsandmorespecicallytheGibbssamplerhavebeenthestandardapproachforBayesianmodels,integratednestedLaplaceapproximations(INLA)methodsofferacomputationallyefcientalternative( Fongetal. ( 2010 )).InthecontextoftheBerksonerrormodelconsideredhere,suchanapproachcouldrepresentaharmonybetweenthedesirablepropertiesofBayesianinferenceandthecomputationalefciencyofthesimplerapproachpresentedhere. 101

PAGE 102

TablesandFigures Table4-1. PoissonResponseSimulationResults:KnownBerksonParameters:Slopeequalto0.02and0.1:Rangeequalto50and400:Scaleequalto70:E(YjX=40)equalto1and48:Average^1,MSE,andvariance^2^1oftheestimatesof1,averageoftheestimatedvariances2^1of^1andcoverage. Method1RangeE(YjX=40)^1AbsoluteBiasMSE^2^1^s2^1Coverage Naive0.025010.0202.50E-063.45E-043.45E-043.25E-040.943PQL0.025010.0203.54E-043.49E-043.49E-043.36E-040.950Naive0.0250480.0201.11E-041.25E-051.25E-056.71E-060.856PQL0.0250480.0201.61E-041.14E-051.14E-051.14E-050.949Naive0.0240010.0195.19E-044.65E-044.64E-044.93E-040.957PQL0.0240010.0203.44E-044.65E-044.65E-044.97E-040.957Naive0.02400480.0202.42E-051.07E-051.07E-051.01E-050.948PQL0.02400480.0201.54E-041.05E-051.05E-051.17E-050.972Naive0.105010.0972.94E-034.38E-044.30E-042.50E-040.862PQL0.105010.0981.81E-033.99E-043.95E-044.15E-040.950Naive0.1050480.0963.57E-032.09E-041.96E-045.17E-060.272PQL0.1050480.1009.27E-059.33E-059.33E-059.57E-050.941Naive0.1040010.0991.02E-036.36E-046.35E-045.70E-040.939PQL0.1040010.1003.86E-046.38E-046.37E-046.23E-040.959Naive0.10400480.1009.80E-056.49E-056.48E-051.16E-050.535PQL0.10400480.1001.35E-044.29E-054.28E-054.01E-050.950 102

PAGE 103

Table4-2. BinomialResponseSimulationResults:KnownBerksonParameters:Slopeequalto0.02and0.1:Rangeequalto50and400:Scaleequalto70:E(YjX=40)equalto0.1and0.5:Average^1,MSE,andvariance^2^1oftheestimatesof1,averageoftheestimatedvariances2^1of^1andcoverage. MethodTrials1RangeE(YjX=40)^1AbsoluteBiasMSE^2^1^s2^1Coverage Naive100.02500.100.0202.71E-043.91E-043.91E-043.65E-040.937PQL100.02500.100.0216.55E-043.97E-043.96E-043.77E-040.945Naive100.02500.500.0203.32E-051.54E-041.54E-041.32E-040.933PQL100.02500.500.0203.77E-041.54E-041.54E-041.40E-040.945Naive100.024000.100.0191.05E-035.13E-045.12E-045.48E-040.955PQL100.024000.100.0199.08E-045.13E-045.13E-045.52E-040.956Naive100.024000.500.0219.20E-041.97E-041.96E-041.95E-040.950PQL100.024000.500.0211.09E-031.98E-041.97E-041.98E-040.952Naive100.10500.100.0954.91E-035.13E-044.89E-043.45E-040.891PQL100.10500.100.0972.72E-035.22E-045.14E-044.98E-040.932Naive100.10500.500.0946.06E-033.33E-042.96E-041.61E-040.824PQL100.10500.500.0981.99E-033.13E-043.09E-042.84E-040.922Naive100.104000.100.1004.23E-046.11E-046.11E-046.49E-040.947PQL100.104000.100.1004.12E-046.13E-046.13E-046.99E-040.964Naive100.104000.500.1001.83E-042.88E-042.88E-042.39E-040.928PQL100.104000.500.1011.05E-032.90E-042.89E-042.80E-040.949Naive500.02500.100.0209.25E-058.11E-058.11E-057.22E-050.929PQL500.02500.100.0202.79E-048.32E-058.31E-057.89E-050.939Naive500.02500.500.0201.74E-043.43E-053.42E-052.64E-050.906PQL500.02500.500.0201.54E-043.42E-053.42E-053.20E-050.946Naive500.024000.100.0204.67E-041.00E-041.00E-041.09E-040.951PQL500.024000.100.0203.08E-041.00E-041.00E-041.11E-040.958Naive500.024000.500.0204.06E-044.02E-054.01E-053.90E-050.945PQL500.024000.500.0215.58E-044.05E-054.02E-054.08E-050.951Naive500.10500.100.0955.00E-032.42E-042.17E-046.86E-050.719PQL500.10500.100.0991.42E-032.02E-042.00E-041.88E-040.929Naive500.10500.500.0946.25E-031.97E-041.58E-043.20E-050.568PQL500.10500.500.1004.89E-041.29E-041.29E-041.31E-040.938Naive500.104000.100.0999.06E-041.57E-041.56E-041.28E-040.918PQL500.104000.100.1001.69E-041.50E-041.50E-041.68E-040.961Naive500.104000.500.1003.76E-049.63E-059.61E-054.75E-050.840PQL500.104000.500.1015.86E-048.79E-058.76E-058.13E-050.951 103

PAGE 104

Table4-3. NaiveInference:Fortheadditionalnormalsimulationswiththerangeequalto200andthescaleequalto105,wereportthemaximumGershgorinradiusandcoveragefornaivecondenceintervalswithknownBerksonparameters. 1RmaxNaiveCoverage 0.12.970.9450.211.860.9370.447.470.9020.6106.800.8670.8189.870.849 Table4-4. NormalResponseSimulationResults:Slopeequalto0.8.Scaleequalto105.Rangeequalto200.Average^1,variance^2^1andMSEoftheestimatesof1,averageoftheestimatedvariances2^1of^1andcoverage. BerksonParametersMethod^1AbsoluteBias^2^1MSE^s2^1Coverage KnownNaive0.7982.20E-033.65E-033.66E-031.89E-030.849KnownMOM0.8054.63E-032.70E-032.72E-032.72E-030.945KnownREML-like0.8054.67E-032.70E-032.72E-032.73E-030.952UnknownNaive0.8055.19E-034.05E-034.08E-031.96E-030.831UnknownMOM0.8003.51E-042.89E-032.89E-032.79E-030.945UnknownREML-Like0.8003.69E-042.88E-032.88E-032.80E-030.951 104

PAGE 105

Figure4-1. Poissonresponsenaivecoverage:ForknownBerksonparameters,aplotofthecoveragefor95%condenceintervalsbasedonthenaivemethodversustheexpectedvalueofYatX=40withaseparatecolorforeachvalueoftherangeparameter.Eachcellcorrespondstoacombinationoftheslopeandscaleparameter. 105

PAGE 106

Figure4-2. PoissonresponsePQL-basedcoverage:ForknownBerksonparameters,aplotofthecoveragefor95%condenceintervalsbasedonthePQLmethodversustheexpectedvalueofYatX=40withaseparatecolorforeachvalueoftherangeparameter.Eachcellcorrespondstoacombinationoftheslopeandscaleparameter. 106

PAGE 107

Figure4-3. Poissonresponsenaive:PQL-basedratio:ForknownBerksonparameters,aplotoftheaverageoftheestimatedvariancebasedonthenaivemethoddividedbytheaverageoftheestimatedvariancebasedonthePQLmethodversustheexpectedvalueofYatX=40withaseparatecolorforeachvalueoftherangeparameter.Eachcellcorrespondstoacombinationoftheslopeandscaleparameter.CurvesareaddedbasedontheapproximaterelationshipderivedinChapter 7 107

PAGE 108

Figure4-4. Binomialresponsenaive:PQL-basedratio:50Trials:ForknownBerksonparameters,aplotoftheaverageoftheestimatedvariancebasedonthenaivemethoddividedbytheaverageoftheestimatedvariancebasedonthePQLmethodversustheexpectedvalueofYatX=40withaseparatecolorforeachvalueoftherangeparameter.Eachcellcorrespondstoacombinationoftheslopeandscaleparameter.CurvesareaddedbasedontheapproximaterelationshipderivedinChapter 7 108

PAGE 109

Figure4-5. PM2.5monitorsinFlorida 109

PAGE 110

CHAPTER5ASYMPTOTICPROPERTIESOFTHEMOMESTIMATORInChapter 3 Section 3.3.1.1 ,forthecasewhen=2Inn,anestimatedgeneralizedsquaresestimatorwasderivedusingamethod-of-momentsestimatorofy.Thederivationoftheestimatorisrestatedherewithadditionaldetailsregardingitsasymptoticproperties.Namely,weprovetheestimatorofyisconsistent,i.e.,convergesinprobabilitytoyand,asaresult,theEGLSestimatorisasymptoticallyequivalenttotheGLSestimator.First,considertheOLSestimatorof, ^OLS=(X0X))]TJ /F5 7.97 Tf 6.58 0 Td[(1X0Y(su).(5)TheOLSestimatorisunbiasedfor; E(^OLS)=(X0X))]TJ /F5 7.97 Tf 6.59 0 Td[(1X0E(Y(su))=(X0X))]TJ /F5 7.97 Tf 6.59 0 Td[(1X0X=. (5) Supposen)]TJ /F5 7.97 Tf 6.58 0 Td[(1(X0X)!V1,andn)]TJ /F5 7.97 Tf 6.59 0 Td[(1X0yX!V2.Asaresult, limn!1Var(^OLS)=limn!1n)]TJ /F5 7.97 Tf 6.59 0 Td[(1(n)]TJ /F5 7.97 Tf 6.59 0 Td[(1X0X))]TJ /F5 7.97 Tf 6.59 0 Td[(1n)]TJ /F5 7.97 Tf 6.59 0 Td[(1X0yX(n)]TJ /F5 7.97 Tf 6.59 0 Td[(1X0X))]TJ /F5 7.97 Tf 6.58 0 Td[(1 (5) =limn!1n)]TJ /F5 7.97 Tf 6.59 0 Td[(1V)]TJ /F5 7.97 Tf 6.58 0 Td[(11V2V)]TJ /F5 7.97 Tf 6.59 0 Td[(11 (5) =0 (5) ) (5) ^OLS=+op(1). (5) 110

PAGE 111

Second,considerthequadraticform, Q1=Y0(I)]TJ /F10 11.955 Tf 11.95 0 Td[(P)Y,(5)whereP=X(X0X))]TJ /F5 7.97 Tf 6.59 0 Td[(1X0.TheexpectedvalueofQ1isE(Q1)=E(Y0(I)]TJ /F10 11.955 Tf 11.95 0 Td[(P)Y)=trace((I)]TJ /F10 11.955 Tf 11.96 0 Td[(P)y)=trace((I)]TJ /F10 11.955 Tf 11.96 0 Td[(P)2I+(I)]TJ /F10 11.955 Tf 11.95 0 Td[(P)(21))=(n)]TJ /F4 11.955 Tf 11.95 0 Td[(p)2+trace((I)]TJ /F10 11.955 Tf 11.95 0 Td[(P)(21)),wherepisthedimensionof.Asaresult, Q2=Y0(I)]TJ /F10 11.955 Tf 11.95 0 Td[(P)Y)]TJ /F4 11.955 Tf 11.95 0 Td[(trace((I)]TJ /F10 11.955 Tf 11.95 0 Td[(P)(21)) n)]TJ /F4 11.955 Tf 11.95 0 Td[(p(5)isunbiasedfor2.Consider^Q2, ^Q2=Y0(I)]TJ /F10 11.955 Tf 11.95 0 Td[(P)Y)]TJ /F4 11.955 Tf 11.95 0 Td[(trace((I)]TJ /F10 11.955 Tf 11.95 0 Td[(P)(^21,OLS)) n)]TJ /F4 11.955 Tf 11.95 0 Td[(p(5)asanestimateofQ2. Theorem5.1. ^Q2isconsistentfor2.Equivalently,wesay^Q2convergesinprobabilityto2or,insymbols, ^Q2p!2,(5)wherep!denotesconvergenceinprobability. Proof. BecauseE(Q2)=2,toproveQ2isconsistentfor2,weprovethevarianceofQ2convergesto0asn!1.Becauseyisasymmetricmatrix,thereexistsadiagonalmatrixandanorthogonalmatrix)]TJ /F1 11.955 Tf 10.27 0 Td[(suchthat y=\012)]TJ /F13 7.97 Tf 23.38 4.93 Td[(0.(5) 111

PAGE 112

Thus, limn!1Var(Q2)=limn!1(n)]TJ /F4 11.955 Tf 11.96 0 Td[(p))]TJ /F5 7.97 Tf 6.59 0 Td[(2trace((I)]TJ /F10 11.955 Tf 11.96 0 Td[(P)y(I)]TJ /F10 11.955 Tf 11.96 0 Td[(P)y) (5) =limn!1(n)]TJ /F4 11.955 Tf 11.96 0 Td[(p))]TJ /F5 7.97 Tf 6.59 0 Td[(2trace((I)]TJ /F10 11.955 Tf 11.96 0 Td[(P)\012)]TJ /F13 7.97 Tf 23.38 4.93 Td[(0(I)]TJ /F10 11.955 Tf 11.95 0 Td[(P)\012)]TJ /F13 7.97 Tf 23.37 4.93 Td[(0) (5) =limn!1(n)]TJ /F4 11.955 Tf 11.96 0 Td[(p))]TJ /F5 7.97 Tf 6.59 0 Td[(2trace()]TJ /F13 7.97 Tf 6.94 4.93 Td[(0(I)]TJ /F10 11.955 Tf 11.95 0 Td[(P)\012)]TJ /F13 7.97 Tf 23.38 4.93 Td[(0(I)]TJ /F10 11.955 Tf 11.96 0 Td[(P)\012) (5) =limn!1(n)]TJ /F4 11.955 Tf 11.96 0 Td[(p))]TJ /F5 7.97 Tf 6.59 0 Td[(2trace((I)]TJ /F10 11.955 Tf 11.96 0 Td[()]TJ /F13 7.97 Tf 6.94 4.94 Td[(0P)]TJ /F6 11.955 Tf 15.35 0 Td[()(I)]TJ /F10 11.955 Tf 11.95 0 Td[()]TJ /F13 7.97 Tf 6.94 4.94 Td[(0P)]TJ /F6 11.955 Tf 15.34 0 Td[()) (5) =limn!1(n)]TJ /F4 11.955 Tf 11.96 0 Td[(p))]TJ /F5 7.97 Tf 6.59 0 Td[(2trace(2)]TJ /F6 11.955 Tf 11.95 0 Td[(2P+PP), (5) whereP=)]TJ /F13 7.97 Tf 6.95 4.34 Td[(0P)]TJ /F1 11.955 Tf 15.34 0 Td[(.Note,trace()]TJ /F13 7.97 Tf 6.94 4.34 Td[(0P)]TJ /F6 11.955 Tf 15.35 0 Td[()=trace(P)=p.BecauseP,P,andaresymmetric, limn!1Var(Q2)=limn!1(n)]TJ /F4 11.955 Tf 11.96 0 Td[(p))]TJ /F5 7.97 Tf 6.59 0 Td[(2trace(2)]TJ /F6 11.955 Tf 11.95 0 Td[(2P+PP), (5) =limn!1(n)]TJ /F4 11.955 Tf 11.96 0 Td[(p))]TJ /F5 7.97 Tf 6.59 0 Td[(2trace(2)]TJ /F6 11.955 Tf 11.95 0 Td[(2P+P2) (5) =limn!1(n)]TJ /F4 11.955 Tf 11.96 0 Td[(p))]TJ /F5 7.97 Tf 6.59 0 Td[(2trace((I+P)2)]TJ /F6 11.955 Tf 11.95 0 Td[(2P) (5) limn!1(n)]TJ /F4 11.955 Tf 11.96 0 Td[(p))]TJ /F5 7.97 Tf 6.59 0 Td[(2((n+p)max1in2ii)]TJ /F6 11.955 Tf 11.95 0 Td[(2pmax1in2ii) (5) =limn!1(n)]TJ /F4 11.955 Tf 11.96 0 Td[(p))]TJ /F5 7.97 Tf 6.59 0 Td[(2((n)]TJ /F4 11.955 Tf 11.96 0 Td[(p)max1in2ii) (5) =0. (5) BecausethevarianceofQ2convergestozero, Q2 2,(5)where denotesconverenceindistribution.Becauseconvergenceindistributiontoaconstantisequivalenttoconvergenceinprobability,theresult Q2p!2(5)istrue.Because^1,OLSp!1, 112

PAGE 113

^21,OLSp!21bythecontinuousmappingtheorem.Again,bythecontinuousmappingtheorem,^Q2)]TJ /F4 11.955 Tf 11.96 0 Td[(Q2p!0.Thus,^Q2p!2bytheconvergingtogethertheorem. BasedontheOLSestimateof1and^Q2,consider ^y=^21,OLS+^Q2Inn(5)asanestimateofy.Because^21,OLSp!21,^21,OLSp!21.BySlutsky'stheorem,^y)]TJ /F6 11.955 Tf 11.96 0 Td[((21+Q2Inn)p!0,and21+Q2Innp!y.Asaresult,^yp!ybytheconvergingtogethertheorem.AnEGLSestimatoristhen ^EGLS=(X0^)]TJ /F5 7.97 Tf 6.59 0 Td[(1yX))]TJ /F5 7.97 Tf 6.59 0 Td[(1X0^)]TJ /F5 7.97 Tf 6.59 0 Td[(1yY(su),(5) 113

PAGE 114

andanestimateofthevarianceis dVar(^EGLS)=(X0^)]TJ /F5 7.97 Tf 6.59 0 Td[(1yX))]TJ /F5 7.97 Tf 6.59 0 Td[(1.(5)Because^yisconsistentfory,^EGLSisasymptoticallyequivalenttotheGLSestimator,^GLS( Maddala ( 1971 )).Becauseyisafunctionof1,theinitialestimatesareusedinaniterativeprocesstoobtainthenalparameterestimates. 1. Update^Q2. 2. Update 5 with^Q2from2and^1,EGLS. 3. Recalculate^EGLSanditsvariance,dVar(^EGLS). 4. Repeatsteps1-3untilchangesintheestimatearesufcientlysmall. 114

PAGE 115

CHAPTER6ASYMPTOTICPROPERTIESOFPSEUDO-LIKELIHOODESTIMATORSINLINEARREGRESSIONMODELS 6.1ModelFormulationSupposetherelationshipbetweenYandXis Y=X+,N(0,(,)) (6) ()YN(X,(,)) (6) whereYisavectorofnobservedvalues,Xisafull(column)ranknqmatrixofcovariates,isanq1vectorofparameters,isap1vectorofparameters,(,)isannnnonsingularpositive-denitevariance-covariancematrix.Thegoalistoestimateanditsstandarderror. 6.2GLSIf(,)=isknown,thenthegeneralizedleastsquaresestimatorof, ^GLS=(X0)]TJ /F5 7.97 Tf 6.58 0 Td[(1X))]TJ /F5 7.97 Tf 6.59 0 Td[(1X0)]TJ /F5 7.97 Tf 6.59 0 Td[(1Y,(6)isthebestlinearunbiasedestimator.Becauseisnormallydistributed, ^GLSN(,(X0)]TJ /F5 7.97 Tf 6.58 0 Td[(1X))]TJ /F5 7.97 Tf 6.59 0 Td[(1)(6)( Khuri ( 2009 )Chapter6).Becausethegeneralizedleastsquaresestimatordependsonthecovariancematrix,anestimatedgeneralizedleastsquares(EGLS)estimatorisfoundusinganestimateof.HerewedescribelikelihoodbasedmethodstoestimateandderivetheasymptoticpropertiesofthesubsequentlydenedEGLSestimator. 115

PAGE 116

6.3EGLSTheestimatorin 6 cannotbeobtainedherebecauseandand,asaresult,(,),arenotknown.However,If~isaconsistentestimatorof,then ~EGLS=(X0~)]TJ /F5 7.97 Tf 6.58 0 Td[(1X))]TJ /F5 7.97 Tf 6.59 0 Td[(1X0~)]TJ /F5 7.97 Tf 6.59 0 Td[(1Y(6)isasymptoticallyequivalenttotheGLSestimatorand ~EGLS N(,(X0)]TJ /F5 7.97 Tf 6.58 0 Td[(1X))]TJ /F5 7.97 Tf 6.59 0 Td[(1).(6)Morespecically,if~and~areconsistentestimatorsofand,then~=(~,~)isaconsistentestimatorof(,)( Maddala ( 1971 ); Magnus ( 1978 )).Thus,inordertoobtainaconsistentestimatorofusingEGLS,consistentestimatorsofandareneeded.Inthenexttwosectionstheseestimatorsarederived. 6.4OLSConsidertheOLSestimateof, ~=(X0X))]TJ /F5 7.97 Tf 6.59 0 Td[(1X0Y.(6)Notice E(~)= (6) and Var(~)=(X0X))]TJ /F5 7.97 Tf 6.59 0 Td[(1X0(,)X(X0X))]TJ /F5 7.97 Tf 6.59 0 Td[(1 (6) =n)]TJ /F5 7.97 Tf 6.59 0 Td[(1(n)]TJ /F5 7.97 Tf 6.58 0 Td[(1X0X))]TJ /F5 7.97 Tf 6.58 0 Td[(1n)]TJ /F5 7.97 Tf 6.58 0 Td[(1X0(,)X(n)]TJ /F5 7.97 Tf 6.58 0 Td[(1X0X))]TJ /F5 7.97 Tf 6.58 0 Td[(1. (6) Assume n)]TJ /F5 7.97 Tf 6.58 0 Td[(1X0X!V1,(6) 116

PAGE 117

whereV1isnonsingular,niteandpositivedeniteandthatV2isanite,positivedenitematrixdenedasn)]TJ /F5 7.97 Tf 6.59 0 Td[(1X0(,)X!V2.Becauseofthenormalityof, n)]TJ /F5 7.97 Tf 6.59 0 Td[(1=2X0 N(0,V2).(6) Proposition6.1. Thefollowingaretrue 1. Var(~)!0 2. ~=+op(1). 3. p n(~)]TJ /F12 11.955 Tf 11.96 0 Td[()N(0,V)]TJ /F5 7.97 Tf 6.58 0 Td[(11V2V)]TJ /F5 7.97 Tf 6.59 0 Td[(11). Proof. limn!1Var(~)=limn!1n)]TJ /F5 7.97 Tf 6.58 0 Td[(1(n)]TJ /F5 7.97 Tf 6.59 0 Td[(1X0X))]TJ /F5 7.97 Tf 6.59 0 Td[(1n)]TJ /F5 7.97 Tf 6.59 0 Td[(1X0(,)X(n)]TJ /F5 7.97 Tf 6.59 0 Td[(1X0X))]TJ /F5 7.97 Tf 6.59 0 Td[(1 (6) =limn!1n)]TJ /F5 7.97 Tf 6.58 0 Td[(1V)]TJ /F5 7.97 Tf 6.59 0 Td[(11V2V)]TJ /F5 7.97 Tf 6.59 0 Td[(11 (6) =0 (6) ) (6) ~=+op(1) (6) Finally, p n(~)]TJ /F12 11.955 Tf 11.96 0 Td[()=p n((X0X))]TJ /F5 7.97 Tf 6.59 0 Td[(1X0Y)]TJ /F12 11.955 Tf 11.95 0 Td[() (6) =(n)]TJ /F5 7.97 Tf 6.58 0 Td[(1X0X))]TJ /F5 7.97 Tf 6.58 0 Td[(1n)]TJ /F5 7.97 Tf 6.59 0 Td[(1=2X0 (6) N(0,,=V)]TJ /F5 7.97 Tf 6.58 0 Td[(11V2V)]TJ /F5 7.97 Tf 6.59 0 Td[(11), (6) because 6 and 6 hold. 117

PAGE 118

6.5REMLConsidertheerrorcontrast (I)]TJ /F10 11.955 Tf 11.95 0 Td[(P)YN(0,(I)]TJ /F10 11.955 Tf 11.95 0 Td[(P)(,)(I)]TJ /F10 11.955 Tf 11.96 0 Td[(P))(6)whereP=X(X0X))]TJ /F5 7.97 Tf 6.59 0 Td[(1X0.Theminustwolog-likelihoodoftheerrorcontrastis )]TJ /F6 11.955 Tf 11.95 0 Td[(2lr=log(j(,)j)+log(jX0(,))]TJ /F5 7.97 Tf 6.59 0 Td[(1Xj)+r0(,))]TJ /F5 7.97 Tf 6.58 0 Td[(1r+c, (6) wherecisaconstantthatdoesnotdependon(,)and r=Y)]TJ /F10 11.955 Tf 11.96 0 Td[(X(X0X))]TJ /F5 7.97 Tf 6.58 0 Td[(1X0Y.(6)( LaMotte ( 2007 ))Supposeweareinterestedinestimatingas ^r()=argmaxlr(,).(6) Assumption1. Suppose (a) Thetruevalues(,)2interior()whereistheparameterspaceof(,)and(,)isidentied. (b) lriscontinuousandtwicedifferentiableon. (c) Theorderofintegrationanddifferentiationofthelog-likelihoodcanbereversedfortherstandsecondderivativesoflrwithrespecttoandthemixedpartialderivativewithrespecttoand. (d) l(,)=)]TJ /F4 11.955 Tf 9.29 0 Td[(n)]TJ /F5 7.97 Tf 6.59 0 Td[(1@2 @0lr(,)j=,=p!I, (6) l(,)=)]TJ /F4 11.955 Tf 9.3 0 Td[(n)]TJ /F5 7.97 Tf 6.59 0 Td[(1@2 @0lr(,)j=,=p!I, (6) I,=I0, (6) whereI,andI,arepositivedeniteandnonsingular. 118

PAGE 119

(e) @3 @i@j@klr(,),1i,j,kp (6) @3 @i@j@klr(,),1iq,1j,kp (6) @3 @i@j@klr(,),1i,jq,1kp (6) aredominatedinaneighborhoodof(,)byanintegrablefunctionofY.Note,theseassumptionsarethestandardassumptionsformaximumlikelihoodestimation.Seeforexample Foutz ( 1977 )or Lehmann&Casella ( 1998 )page463. Proposition6.2. Ifisknownthentherestrictedmaximumlikelihoodestimate^r()satises ^r()=+op(1) (6) ^r()=+Op(n)]TJ /F5 7.97 Tf 6.59 0 Td[(1=2) (6) Beforeprovingtheproposition,werestateTheorem3.1from Cressie&Lahiri ( 1993 )asthisresultwillbeusedtoprovetheproposition. Theorem6.1. Assume 1. (,)iscontinuousandtwicedifferentiableon. 2. ThereexistnonrandommatricesBnsuchthatjjB)]TJ /F5 7.97 Tf 6.59 0 Td[(1njjconvergesuniformlyto0and B)]TJ /F5 7.97 Tf 6.59 0 Td[(1n@2 @0lr(,)j=,=(B)]TJ /F5 7.97 Tf 6.58 0 Td[(1n)0p!I,,(6) 3. Forallc>0,>0, (a) supjjB)]TJ /F5 7.97 Tf 6.59 0 Td[(1nBn)]TJ /F10 11.955 Tf 11.95 0 Td[(Ijj:jj(Bn)0()]TJ /F12 11.955 Tf 11.95 0 Td[(0)jjc;,02convergesuniformlyto0;and (b) LetM=(01,...,0k),0i2,whereistheportionoftheparameterspacewherelies,andA(M)bethematrixwith(i,j)thelementequalto@2 @ijlr(,)j=0i,=, 119

PAGE 120

1i,jp. PsupjjB)]TJ /F5 7.97 Tf 6.58 0 Td[(1n(A(M))]TJ /F24 10.909 Tf 17.51 7.38 Td[(@2 @0lr(,)j=,=)(B)]TJ /F5 7.97 Tf 6.59 0 Td[(1n)0jj:jjB0n()]TJ /F40 10.909 Tf 10.91 0 Td[(0i)jj(6)convergesuniformlyto0.Let^nbetheREMLestimatorbasedontherstnobservations.Then,(Bn))]TJ /F5 7.97 Tf 6.59 0 Td[(1(^n)]TJ /F12 11.955 Tf 11.96 0 Td[() N(0,I,).Usingthistheorem,weprovetheproposition. ProofofProposition 6.2 LetBn=n1=2Ipp.TheassumptionsoutlinedinAssumption 1 implytheassumptionsinTheorem 6.1 .Asaresult, p n(^r())]TJ /F12 11.955 Tf 11.95 0 Td[() N(0,I)]TJ /F5 7.97 Tf 6.58 0 Td[(1,).(6)Thus, ^r()=+op(1),and (6) ^r()=+Op(n)]TJ /F5 7.97 Tf 6.58 0 Td[(1=2). (6) 6.6PseudoRestrictedLikelihoodLetthepseudolikelihoodfunction( Gong&Samaniego ( 1981 ); Parke ( 1986 ))bedenedbypluggingtheOLSestimateofintothelikelihoodequationfortheerrorcontrast; )]TJ /F6 11.955 Tf 11.95 0 Td[(2lr(,~)=log(j(,~)j)+log(jX0(,~))]TJ /F5 7.97 Tf 6.58 0 Td[(1Xj)+r0(,~))]TJ /F5 7.97 Tf 6.59 0 Td[(1r+c. (6) Considertheestimateofobtainedbymaximizingthepseudolikelihood, ^r(~)=argmaxlr(,~).(6) Theorem6.2. UndertheassumptionsinAssumption 1 120

PAGE 121

1. p n(^r(~))]TJ /F12 11.955 Tf 11.37 0 Td[() N(0,I)]TJ /F5 7.97 Tf 6.59 0 Td[(1,+I)]TJ /F5 7.97 Tf 6.58 0 Td[(1,I,,I,I)]TJ /F5 7.97 Tf 6.59 0 Td[(1,),withthematricesinthecovariancedenedinAssumption 1 2. ^r(~)p!. Proof. Let l(,)=n)]TJ /F5 7.97 Tf 6.59 0 Td[(1@ @lr(,)j=,=.(6)Undertheassumedmodel,forxed,theequationl(,)=0hasauniqueconsistentroot^r().Thatis, l(^r(),)=0.(6)Expandingl(^r(),)aroundthetruevalueyields 0=l(^r(),)=l(,)+@ @l(,)(^r())]TJ /F12 11.955 Tf 11.95 0 Td[()+op(n)]TJ /F5 7.97 Tf 6.59 0 Td[(1=2) (6) ))]TJ /F6 11.955 Tf 8.08 2.66 Td[(l(,)=@ @l(,)(^r())]TJ /F12 11.955 Tf 11.95 0 Td[()+op(n)]TJ /F5 7.97 Tf 6.59 0 Td[(1=2). (6) Byassumption, l(,)=)]TJ /F10 11.955 Tf 9.3 0 Td[(I,+op(1). (6) Asaresult, )]TJ /F6 11.955 Tf 10.73 2.66 Td[(l(,)=()]TJ /F10 11.955 Tf 9.3 0 Td[(I,+op(1))(^r())]TJ /F12 11.955 Tf 11.95 0 Td[()+op(n)]TJ /F5 7.97 Tf 6.58 0 Td[(1=2) (6) =)]TJ /F10 11.955 Tf 9.3 0 Td[(I,(^r())]TJ /F12 11.955 Tf 11.96 0 Td[()+op(1)(^r())]TJ /F12 11.955 Tf 11.96 0 Td[()+op(n)]TJ /F5 7.97 Tf 6.58 0 Td[(1=2) (6) )]TJ 9.3 9.44 Td[(p nl(,)=)]TJ 9.3 9.44 Td[(p nI,(^r())]TJ /F12 11.955 Tf 11.95 0 Td[()+op(1)p n(^r())]TJ /F12 11.955 Tf 11.95 0 Td[()+op(1) (6) =)]TJ 9.3 9.44 Td[(p nI,(^r())]TJ /F12 11.955 Tf 11.95 0 Td[()+op(1). (6) 121

PAGE 122

Byproposition 6.2 op(1)p n(^r())]TJ /F12 11.955 Tf 11.95 0 Td[()=op(1)Op(1)=op(1) (6) )p n(^r())]TJ /F12 11.955 Tf 11.95 0 Td[()=p nI)]TJ /F5 7.97 Tf 6.58 0 Td[(1,l(,)+op(1) (6) Bydifferentiating 6 withrespecttoandevaluatingatthetruevalue, p n@ @^r()j==p nI)]TJ /F5 7.97 Tf 6.58 0 Td[(1,l(,)+op(1). (6) TheasymptoticdistributionfollowsfromtheTaylorseriesapproximation: p n(^r(~))]TJ /F40 10.909 Tf 10.91 0 Td[()=p n(^r())]TJ /F40 10.909 Tf 10.9 0 Td[()+p n@ @^r()(~)]TJ /F40 10.909 Tf 10.91 0 Td[()+op(1) (6) =p n(^r())]TJ /F40 10.909 Tf 10.9 0 Td[()+p nI)]TJ /F5 7.97 Tf 6.59 0 Td[(1,l(,)(~)]TJ /F40 10.909 Tf 10.91 0 Td[()+op(1)(~)]TJ /F40 10.909 Tf 10.91 0 Td[()+op(1) (6) =p n(^r())]TJ /F40 10.909 Tf 10.9 0 Td[()+p nI)]TJ /F5 7.97 Tf 6.59 0 Td[(1,l(,)(~)]TJ /F40 10.909 Tf 10.91 0 Td[()+op(1). (6) Byassumption,l=I,+op(1).Thus, p n(^r(~))]TJ /F12 11.955 Tf 11.96 0 Td[()=p n(^r())]TJ /F12 11.955 Tf 11.96 0 Td[())]TJ 11.95 9.44 Td[(p nI)]TJ /F5 7.97 Tf 6.58 0 Td[(1,I,(~)]TJ /F12 11.955 Tf 11.96 0 Td[()+p nI)]TJ /F5 7.97 Tf 6.59 0 Td[(1,op(1)(~)]TJ /F12 11.955 Tf 11.96 0 Td[()+op(1). (6) Byproposition 6.1 ,p n(~)]TJ /F12 11.955 Tf 11.96 0 Td[()isOp(1).Asaresult, I)]TJ /F5 7.97 Tf 6.59 0 Td[(1,op(1)p n(~)]TJ /F12 11.955 Tf 11.96 0 Td[()=I)]TJ /F5 7.97 Tf 6.58 0 Td[(1,op(1)Op(1)=I)]TJ /F5 7.97 Tf 6.58 0 Td[(1,op(1)=op(1),(6)and p n(^r(~))]TJ /F12 11.955 Tf 11.96 0 Td[()=p n(^r())]TJ /F12 11.955 Tf 11.96 0 Td[())]TJ 11.95 9.44 Td[(p nI)]TJ /F5 7.97 Tf 6.58 0 Td[(1,I,(~)]TJ /F12 11.955 Tf 11.96 0 Td[()+op(1). (6) Asamaximumlikelihoodestimator,^r()isasymptoticallyefcient.Inaddition,p n(~)]TJ /F12 11.955 Tf 12.45 0 Td[()hasasymptoticmeanzero.Asaresult,p n(^())]TJ /F12 11.955 Tf 12.45 0 Td[()andp n(~)]TJ /F12 11.955 Tf 12.44 0 Td[()are 122

PAGE 123

asymptoticallyindependent( Parke ( 1986 )).Thus, p n(^r(~))]TJ /F12 11.955 Tf 11.95 0 Td[() N(0,I)]TJ /F5 7.97 Tf 6.59 0 Td[(1,+I)]TJ /F5 7.97 Tf 6.59 0 Td[(1,I,,I,I)]TJ /F5 7.97 Tf 6.59 0 Td[(1,)(6)Finally,because )]TJ /F10 11.955 Tf 11.96 0 Td[(I)]TJ /F5 7.97 Tf 6.59 0 Td[(1,I,(~)]TJ /F12 11.955 Tf 11.95 0 Td[()=op(1), (6) thefollowingistrueandtheresultisproven; ^r(~)=^r())]TJ /F10 11.955 Tf 11.95 0 Td[(I)]TJ /F5 7.97 Tf 6.59 0 Td[(1,I,(~)]TJ /F12 11.955 Tf 11.96 0 Td[()+op(n)]TJ /F5 7.97 Tf 6.59 0 Td[(1=2) (6) =+op(1))]TJ /F10 11.955 Tf 11.96 0 Td[(I)]TJ /F5 7.97 Tf 6.59 0 Td[(1,I,(~)]TJ /F12 11.955 Tf 11.95 0 Td[()+op(n)]TJ /F5 7.97 Tf 6.59 0 Td[(1=2) (6) =+op(1). (6) 6.7EGLSEstimatorAsaresultofTheorem 6.2 ,^=(^r(~),~)isaconsistentestimatorof.Thus, ^EGLS=(X0^)]TJ /F5 7.97 Tf 6.58 0 Td[(1X))]TJ /F5 7.97 Tf 6.59 0 Td[(1X0^)]TJ /F5 7.97 Tf 6.59 0 Td[(1Y(6)isasymptoticallyequivalenttotheGLSestimatorand ^EGLS N(,(X0)]TJ /F5 7.97 Tf 6.58 0 Td[(1X))]TJ /F5 7.97 Tf 6.59 0 Td[(1).(6) 123

PAGE 124

CHAPTER7FUTUREWORK 7.1Spatio-TemporalMisalignmentTheinducedBerksonerrormodelholdsinotherapplicationswhenmissingcovariatesarepredicted.Forexample,ifspatio-temporalkrigingisusedtopredictavariablemisalignedinspaceandtime,thentheproposedmethodologycanbeusedtoaccountfortheadditionaluncertainty.Regardlessofthenatureofspatio-temporalcovariancestructure,e.g.,separableversusnotseparableinspaceandtime,aslongastheBerksonerrorrelationshipholdsbetweentheunobservedcovariateandconditionalmean,theproposedmethodologywillstillapply.Futureworkexiststoassesstheperformanceoftheproposedmethodologyinpractice. 7.2ClassicalMeasurementErrorInChapter 4 ,theproposedmethodologiesweredevelopedforthegeneralBerksonerrorproblem.Forclassicalmeasurementerror(CME),regressioncalibration( Carrolletal. ( 2006 ),Chapter4)isusedtoinduceaBerksonerrormodel.However,uncertaintyestimatesarecalculatedusingabootstrap( Carrolletal. ( 2006 ),Chapter4)oranasymptoticsandwichestimator( Carrolletal. ( 2006 ),AppendixB).WeproposeusingthemethodsinChapter 4 forCMEproblems.ConsiderthefollowingargumentforaPoissonregressionmodel. 124

PAGE 125

7.2.1PoissonRegressionwithCMEConsiderthemeasurementerrormodel, YPoisson() (7) log()=01n1+1X (7) X=x+x (7) W=X+w (7) wN(0,2I) (7) xN(0,2xI). (7) SupposeYandWareobservedandx,2,and2xareknown.ConsidertheconditionaldistributionofXjW: XjWN(XjW=x+2x 2x+2(W)]TJ /F12 11.955 Tf 11.96 0 Td[(x),2RCI),(7)where 2RC=2x)]TJ /F7 11.955 Tf 26.47 8.09 Td[(4x 2x+2.(7)Thatis, XjW=XjW+v,vN(0,2RCI).(7)Considertheinducedmodel, YPoisson() (7) log()=01+1XjW+1v (7) vN(0,2RCI). (7) NoticethemodelnowhastheBerksonerrormodelformasdenedinChapter 4 .Thus,themethodsdescribedinChapter 4 canbeusedinthecaseofclassicalmeasurementerror. 125

PAGE 126

7.2.2SimulationLetn=100,x=40,2=4,2x=4,1=0.2and0=log(32))]TJ /F7 11.955 Tf 12.69 0 Td[(140.Wesimulated3000datasetsfromthemodelin 7 andusedthePQLmethodandthenaivemethodtoconductinferenceabouttheparameter1.Theresults,Table 7-1 ,show,forthePQLmethod,unbiasedestimatesoftheregressionparameteranditsuncertaintyandcoverageisnear95%.Asexpected,forthenaivemethod,theuncertaintyisunderestimatedandcoverageisnear68%. 7.3BerksonErrorInadditiontousingtheproposedmethodologyinChapter 4 toconductinferenceinCMEproblems,weproposethemethodscanbeusedinconjunctionwithmaximumrestrictedpseudo-likelihoodtoconductinferenceinBerksonerrorproblemswithunknownBerksonparameters.Inparticular,weconsiderthecasewhentheBerksonparametersdenethecovarianceofthemeasurementerrorparameters.ConsiderthefollowingargumentforthePoissonregressionmodel. 7.3.1PoissonRegression:UnknownBerksonParameterConsiderthemodel YPoisson() (7) log()=01n1+1X (7) X=W+x (7) xN(0,2xI). (7) SupposeWisknownand2xisnot.BasedonthePQLtechniquesdescribedin Breslow&Clayton ( 1993 ),maximumrestrictedpseudo-likelihoodcanbeusedtosimultaneouslyconductinferenceabout1andestimate2x. 126

PAGE 127

7.3.2SimulationLetn=100.FixW=4011001+v,vN(0,I).Let2x=4,1=0.2and0=log(32))]TJ /F7 11.955 Tf 12.63 0 Td[(140.Wesimulated1000datasetsfromthemodelin 7 andusedthePQLmethodandthenaivemethodtoconductinferenceabouttheparameter1.Inaddition,weusedmaximumrestrictedpseudo-likelihoodtoestimate2x.Theresults,Table 7-2 ,show,forthePQLestimate,unbiasedestimatesof1andit'suncertaintyandcoveragenear95%.Again,asexpected,naivemethodsunderestimatethevariabilityandcoverageisnear53%.ThemedianoftheestimatedBerksonerrorvariance2x=4isequalto4.09. 127

PAGE 128

7.4TablesandFigures Table7-1. Classicalmeasurementerrorsimulationresults:Fortheclassicalmeasurementerrorsimulation,wereportthemeanoftheestimated1values,^1,thevarianceoftheestimated1values,^2^1,themeanoftheestimatedvariancesof^1,^s2^1,andcoveragefor95%forthenaiveandPQLmethod. Method^1^2^1^s2^1Coverage PQL0.1980.0005800.0005860.94Naive0.2000.0006310.0001500.68 Table7-2. BerksonErrorsimulationresults:FortheBerksonerrorsimulation,wereportthemeanoftheestimated1values,^1,thevarianceoftheestimated1values,^2^1,themeanoftheestimatedvariancesof^1,^s2^1,andcoveragefor95%forthenaiveandPQLmethod.Inaddition,themedianoftheestimatedvarianceparameter^2xisreported. Method^1^2^1^s2^1Median^2xCoverage PQL0.1952.25e-32.27e-34.090.95Naive0.1962.40e-33.40e-4NA0.53 128

PAGE 129

APPENDIX:SUPPLEMENTARYMATERIALS A.1ShowingQ2=Q3Q3=X0o0Muo)]TJ /F5 7.97 Tf 6.59 0 Td[(1ooXo)]TJ /F7 11.955 Tf 11.95 0 Td[(X0o0Muo)]TJ /F5 7.97 Tf 6.59 0 Td[(1oo1n1 (A)=X0o0M[uo)]TJ /F5 7.97 Tf 6.59 0 Td[(1ooXo)]TJ /F10 11.955 Tf 11.95 0 Td[(uo)]TJ /F5 7.97 Tf 6.59 0 Td[(1oo1n1] (A)=X0o0M[uo)]TJ /F5 7.97 Tf 6.59 0 Td[(1oo(Xo)]TJ /F10 11.955 Tf 11.95 0 Td[(1n1)] (A)Q2=X0o0MXo (A)=X0o0M[1m1+uo)]TJ /F5 7.97 Tf 6.59 0 Td[(1oo(Xo)]TJ /F10 11.955 Tf 11.96 0 Td[(1n1)] (A)=X0o0M[uo)]TJ /F5 7.97 Tf 6.59 0 Td[(1oo(Xo)]TJ /F10 11.955 Tf 11.95 0 Td[(1n1)] (A)ThelastequalityisduetothefactthatM1m1=0.Thus,youcanseeQ2=Q3. A.2EigenvaluesofInducedCovarianceMatricesSupposethematrixnnmatrixcanbewrittenas =21x+2yI,(A)wherexisapositive-denitecovariancematrix. TheoremA.1. Alleigenvaluesofarecontainedintheinterval[2y,2y+Rmax],whereRi=Pi6=jji,jj,i=1,...,n,andRmax=maxiRi. Proof. Becauseisapostiivedenitecovariancematrix,alleigenvaluesofarerealandpositive.BytheGershgorinDiskTheorem,alleigenvaluesarecontainedinatleastoneinterval[i,i)]TJ /F4 11.955 Tf 11.43 0 Td[(Ri,i,i+Ri],i=1,...,n,wherei,jisthe(i,j)thelementofthematrix(seeTheorem7.2.1 Golub&vanLoan ( 1996 ),p.320). 129

PAGE 130

Becausealldiagonalelementsofxarepositive,alldiagonalelementsofaregreaterthanorequalto2y.Thus,alleigenvaluesarecontainedinatleastoneinterval[2y)]TJ /F4 11.955 Tf 12.35 0 Td[(Ri,2y+Ri],i=1,...,n.Inotherwords,alleigenvalueslieintheinterval[2y)]TJ /F4 11.955 Tf 11.96 0 Td[(Rmax,2y+Rmax].Asaresult,toprovetheresult,weprovethatalleigenvaluesofaregreaterthanorequalto2y.ByWeyl'sinequality( Bhatia ( 1997 ),sectionIII.2,p.62),because,21xand2yIareallpositive-deniteHermitianmatrices,wecanwrite 2y+min(21x)min()2y+max(21x)2ymin()2y+max(21x)2ymin() (A) wheremin(21x)andmax(21x)arethesmallestandlargesteigenvaluesofthematrix21xandmin()isthesmallesteigenvalueof.Thus, A andthefactthatalleigenvalueslieintheinterval[2y)]TJ /F4 11.955 Tf 12.11 0 Td[(Rmax,2y+Rmax]togetherimplyalleigenvaluesoflieintheinterval[2y,2y+Rmax]. A.3ExpectedHealthCountsToaccountfordemographicvariation,theexpectednumberofhealtheventswerecalculatedforeachcounty.Floridawaspartitionedinto80mutuallyexclusivegroupsbytakingallcombinationsof 10CategoriesofAge(inYears):f[0,4],[5,12],[13,17],[18,24],[25,34],[35,44],[45,54],[55,64],[65,74],[75,1)g Gender:fMale,Femaleg Race:fHispanic,Black,White,OthergLeti=1,...,80indexthe80groups.Letp(s)i,i=1,...,80bethepopulationsizesofthe80groupsinthestandardpopulation,here,theentirestateofFlorida. 130

PAGE 131

Foreachgroup,lete(s)i,i=1,...,80bethenumberofcasesinthestandardpopulationduringthetimeperiodofinterest.Thecruderatewascalculatedforeachofthe80groupsas (s)i=e(s)i p(s)i(A)Foreachofthe67countiesinFlorida,letthenumberofpeopleandnumbereventsineachofthegroupsbe eci,i=1,...,80,c=1,...,67(A)and pci,i=1,...,80,c=1,...,67,(A)respectively.LetY(c)=P80i=1eciandE(c)=P80i=1(s)ipcidenotetheobservedandexpectednumberofcasesincountyc,respectively,c=1,...,67. 131

PAGE 132

FigureA-1. Simulationlocations 132

PAGE 133

FigureA-2. Correlationfunctions 133

PAGE 134

FigureA-3. Poissonmeans 134

PAGE 135

FigureA-4. Binomialmeans 135

PAGE 136

FigureA-5. Poissonresponsenaivevariance:ForknownBerksonparameters,aplotoftheaverageoftheestimatedvariancebasedonthenaivemethoddividedbytheestimatedvariabilityoftheslopeparameterversustheexpectedvalueofYatX=40withaseparatecolorforeachvalueoftherangeparameter.Eachcellcorrespondstoacombinationoftheslopeandscaleparameter. 136

PAGE 137

FigureA-6. PoissonresponsePQL-basedvariance:ForknownBerksonparameters,aplotoftheaverageoftheestimatedvariancebasedonthePQLmethoddividedbytheestimatedvariabilityoftheslopeparameterversustheexpectedvalueofYatX=40withaseparatecolorforeachvalueoftherangeparameter.Eachcellcorrespondstoacombinationoftheslopeandscaleparameter. 137

PAGE 138

FigureA-7. Binomialresponsenaivevariance:50Trials:ForknownBerksonparameters,aplotoftheaverageoftheestimatedvariancebasedonthenaivemethoddividedbytheestimatedvariabilityoftheslopeparameterversustheexpectedvalueofYatX=40withaseparatecolorforeachvalueoftherangeparameter.Eachcellcorrespondstoacombinationoftheslopeandscaleparameter. 138

PAGE 139

FigureA-8. Binomialresponsenaivecoverage:50Trials:ForknownBerksonparameters,aplotofthecoveragefor95%condenceintervalsbasedonthenaivemethodversustheexpectedvalueofYatX=40withaseparatecolorforeachvalueoftherangeparameter.Eachcellcorrespondstoacombinationoftheslopeandscaleparameter. 139

PAGE 140

FigureA-9. BinomialresponsePQL-basedvariance:50Trials:ForknownBerksonparameters,aplotoftheaverageoftheestimatedvariancebasedonthePQLmethoddividedbytheestimatedvariabilityoftheslopeparameterversustheexpectedvalueofYatX=40withaseparatecolorforeachvalueoftherangeparameter.Eachcellcorrespondstoacombinationoftheslopeandscaleparameter. 140

PAGE 141

FigureA-10. BinomialresponsePQL-basedcoverage:50Trials:ForknownBerksonparameters,aplotofthecoveragefor95%condenceintervalsbasedonthePQLmethodversustheexpectedvalueofYatX=40withaseparatecolorforeachvalueoftherangeparameter.Eachcellcorrespondstoacombinationoftheslopeandscaleparameter. 141

PAGE 142

FigureA-11. NASAsatelliteimagefromApril29th,2007 FigureA-12. NASAsatelliteimagefromApril30th,2007 142

PAGE 143

FigureA-13. FinegridlaidoverFlorida 143

PAGE 144

REFERENCES Agresti,A.(2002).CategoricalDataAnalysis.NewYork,NewYork:JohnWileyandSons,2nded. Banerjee,S.,&Gelfand,A.(2002).Prediction,interpolationandregressionforspatiallymisaligneddata.TheIndianJournalofStatistics,SeriesA,64(2),227. Bhatia,R.(1997).Matrixanalysis,vol.169.NewYork,NewYork:SpringerVerlag. Breslow,N.,&Clayton,D.(1993).Approximateinferenceingeneralizedlinearmixedmodels.JournaloftheAmericanStatisticalAssociation,88(421),9. Buonaccorsi,J.(2009).Measurementerror:models,methodsandapplications.BocaRaton,FL:Chapman&Hall/CRC. Carroll,R.,Ruppert,D.,Stefanski,L.,&Crainiceanu,C.(2006).MeasurementErrorinNonlinearModels:AModernPerspective.BocaRaton,FL:Chapman&Hall/CRC,2nded. Cressie,N.,&Lahiri,S.(1993).Theasymptoticdistributionofremlestimators.Journalofmultivariateanalysis,45(2),217. FloridaDepartmentofEnvironmentalProtection(2007).Healthadvisory:May10th,2007.Http://www.dep.state..us/mainpage/em/2007/wildre/health/Advisory 051007.pdf. Fong,Y.,Rue,H.,&Wakeeld,J.(2010).Bayesianinferenceforgeneralizedlinearmixedmodels.Biostatistics,11(3),397. Foutz,R.(1977).Ontheuniqueconsistentsolutiontothelikelihoodequations.JournaloftheAmericanStatisticalAssociation,72(357),147. Gelfand,A.,Zhu,L.,&Carlin,B.(2001).Onthechangeofsupportproblemforspatio-temporaldata.Biostatistics,2(1),31. GeorgiaForestryCommision(2007).Thehistoric2007georgiawildres:Learningfromthepast-planningforthefuture.Http://www.wildrelessons.net/documents/Historic 2007 GA Wildres.pdf. Golub,G.,&vanLoan,C.(1996).MatrixComputation.Baltimore,MD:JohnsHopkinsUniversityPress,3ed. Gong,G.,&Samaniego,F.(1981).Pseudomaximumlikelihoodestimation:theoryandapplications.TheAnnalsofStatistics,9(4),861. Gryparis,A.,Paciorek,C.,Zeka,A.,Schwartz,J.,&Coull,B.(2009).Measurementerrorcausedbyspatialmisalignmentinenvironmentalepidemiology.Biostatistics,10(2),258. 144

PAGE 145

Khuri,A.(2009).Linearmodelmethodology.BocaRaton,FL:Chapman&Hall/CRC. LaMotte,L.(2007).Adirectderivationoftheremllikelihoodfunction.StatisticalPapers,48(2),321. Leem,J.-H.,Kaplan,B.M.,Shim,Y.K.,Pohl,H.R.,Gotway,C.A.,Bullard,S.M.,Rogers,J.F.,Smith,M.M.,&Tylenda,C.A.(2006).Exposurestoairpollutantsduringpregnancyandpretermdelivery.EnvironmentalHealthPerspectives,114,905. Lehmann,E.,&Casella,G.(1998).Theoryofpointestimation.NewYork:SpringerVerlag,2nded. Lopiano,K.,Young,L.,&Gotway,C.(2011).Acomparisonoferrorsinvariablesmethodsforuseinregressionmodelswithspatiallymisaligneddata.StatisticalMethodsinMedicalResearch,20(1),29. Maddala,G.(1971).Generalizedleastsquareswithanestimatedvariancecovariancematrix.Econometrica:JournaloftheEconometricSociety,39(1),23. Madsen,L.,Ruppert,D.,&Altman,N.(2008).Regressionwithspatiallymisaligneddata.Environmetrics,19(5),453. Madsen,L.J.(2004).RegressionwithSpatiallyMisalignedData.Ph.D.thesis,CornellUniversity,Ithaca,NY. Magnus,J.(1978).Maximumlikelihoodestimationoftheglsmodelwithunknownparametersinthedisturbancecovariancematrix.Journalofeconometrics,7(3),281. McCullagh,P.,&Nelder,J.(1989).MonographsonStatisticsandAppliedProbability37:GeneralizedLinearModels.BocaRaton,FL:Chapman&Hall/CRC,2nded. Mugglin,A.,&Carlin,B.(1998).Hierarchicalmodelingingeographicinformationsystems:Populationinterpolationoverincompatiblezones.JournalofAgricultural,BiologicalandEnvironmentalStatistics,3(2),111. Mugglin,A.,Carlin,B.,&Gelfand,A.(2000).Fullymodel-basedapproachesforspatiallymisaligneddata.JournaloftheAmericanStatisticalAssociation,95(451),877. Mugglin,A.,Carlin,B.,Zhu,L.,&Conlon,E.(1999).Bayesianarealinterpolation,estimationandsmoothing:Aninferentialapproachforgeographicinformationsystems.EnvironmentandPlanning,Ser.A,31(8),1337. NationalCenterforHealthStatistics(1998).Classicationofdiseasesandinjuries.Ftp://ftp.cdc.gov/pub/Health Statistics/NCHS/Publications/ICD-9/ucod.txt. Paciorek,C.,Yanosky,J.,Puett,R.,Laden,F.,&Suh,H.(2009).Practicallarge-scalespatio-temporalmodelingofparticulatematterconcentrations.TheAnnalsofAppliedStatistics,3(1),370. 145

PAGE 146

Parke,W.(1986).Pseudomaximumlikelihoodestimation:Theasymptoticdistribution.TheAnnalsofStatistics,14(1),355. Rogers,J.F.,Thompson,S.J.,Addy,C.L.,McKeown,R.E.,Cowen,D.J.,&DeCoulfe,P.(2000).Theassociationofverylowbirthweightwithexposurestoenvironmentalsulfurdioxideandtotalsuspendedparticulates.AmericanJournalofSocialEpidemiol-ogy,151(6),602. Samet,J.,Dominici,F.,Curriero,F.,Coursac,I.,&Zeger,S.(2000).Fineparticulateairpollutionandmortalityin20u.s.cities,1987-1994.TheNewEnglandJournalofMedicine,343,1742. Schabenberger,O.,&Gotway,C.(2005).StatisticalMethodsforSpatialDataAnalysis.BocaRaton,FL:Chapman&Hall/CRC,1sted. Schouten,J.,Vonk,J.,&deGraaf,A.(1996).Shorttermeffectsofairpollutiononemergencyhospitaladmissionsforrespiratorydisease:resultsoftheapheaprojectintwomajorcitiesinthenetherlands,1977-89.JournalofEpidemiologicalCommunityHealth,50:Suppl.1,s22. Szpiro,A.,Sheppard,L.,&Lumley,T.(2009).Efcientmeasurementerrorcorrectionwithspatiallymisaligneddata.UniversityofWashingtonBiostatisticsWorkingPaperSeries,Paper350. Szpiro,A.,Sheppard,L.,&Lumley,T.(2011).Efcientmeasurementerrorcorrectionwithspatiallymisaligneddata.Biostatistics,12(4),610. TheStateofFloridaStateEmergencyResponseTeam(2007).Statewidedisasterresponsemorningbrieng:May9th,2007.Http://www.dep.state..us/mainpage/em/2007/wildre/response/20070509 brieng.pdf. Thurston,S.,Spiegelman,D.,&Ruppert,D.(2003).Equivalenceofregressioncalibrationmethodsinmainstudy/externalvalidationstudydesigns.Journalofstatisticalplanningandinference,113(2),527. U.S.EnvironmentalProtectionAgency(1999a).Environmentalmonitoringandassessmentprogram:Mid-atlanticstreamschemistrydata.Http://www.epa.gov/emap/html/data/surfwatr/data/mastreams/9396/wchem/chmval.txt. U.S.EnvironmentalProtectionAgency(1999b).Environmentalmonitoringandassessmentprogram:Mid-atlanticstreamswatershedcharacteristics.Http:/www.epa.gov/emap/html/data/surfwatr/data/mastreams/9396/location/watchr.txt. U.S.EnvironmentalProtectionAgency(2007).Airqualitysystem(aqs)datafordownloading.Http://www.epa.gov/ttn/airs/airsaqs/detaildata/downloadaqsdata.htm. Waller,L.A.,&Gotway,C.A.(2004).AppliedSpatialStatisticsforPublicHealthData.NewYork:JohnWileyandSons. 146

PAGE 147

Woodward,M.(2004).Epidemiology:StudyDesignandDataAnalysis.BocaRaton,Florida:ChapmanandHall/CRC,2nded. Young,L.,&Gotway,C.(2007).Linkingspatialdatafromdifferentsources:Theeffectofchangeofsupport.StochasticEnvironmentalResearchandRiskAssessment,21,589. Young,L.,Gotway,C.,Kearney,G.,&DuClos,C.(2009a).Assessinguncertaintyinsupport-adjustedspatialmisalignmentproblems.CommunicationsinStatisticsTheoryandMethods,38(16-17),3249. Young,L.,Gotway,C.,Yang,J.,Kearney,G.,&DuClos,C.(2009b).Linkinghealthandenvironmentaldataingeographicalanalysis:Itssomuchmorethancentroids.SpatialandSpatio-temporalEpidemiology,1(1),73. Zhang,C.,Tang,Y.,Xu,X.,&Kiely,G.(2011).Towardsspatialgeochemicalmodelling:useofgeographicallyweightedregressionformappingsoilorganiccarboncontentsinireland.AppliedGeochemistry,26(7),1239. Zhu,L.,Carlin,B.,&Gelfand,A.(2003).Hierarchicalregressionwithmisalignedspatialdata:relatingambientozoneandpediatricasthmaervisitsinatlanta.Environmetrics,14(5),537. 147

PAGE 148

BIOGRAPHICALSKETCH KennethKennyLopianowasbornin1986.KennyreceivedaB.S.inMathematicsandaB.S.inStatisticsfromtheUniversityofFloridainMay2007.InMay2009,KennyreceivedhisMastersofStatisticsfromtheUniversityofFlorida.UponcompletionofhisPh.D.programinthesummerof2012,KennybecameapostdoctoralfellowattheStatisticsandAppliedMathematicalSciencesInstitute(SAMSI)inResearchTrianglePark,NC. 148