chronic exposure to 
 ambient pm2.5 is linked to increased 
 morbidity and mortality in many epidemiological studies and results in high population burden , making it a 
 large public health concern . hence quantifying accurate air pollution 
 exposure has become paramount and has prompted different approaches 
 to estimate chronic pm2.5 levels across space and time . as our awareness of the impact of air pollution has increased , 
 so has the interdisciplinary nature of exposure assessment . researchers 
 from these disciplines range from air pollution scientists to epidemiologists 
 to risk assessors who are all involved in better understanding air 
 pollution processes and its health effects . disciplines also extend 
 to cost - benefit analysts , policy makers and regulators whose goals 
 are air pollution abatement through policy to efficiently diminish 
 its burden on the population . because of the wide range of groups 
 involved there is a critical need for methods that are accurate in 
 estimating chronic levels of pm2.5 and are both accessible 
 and interpretable by a wide audience . it is this wide audience which 
 we are keeping in mind in advancing methods used to estimate chronic 
 pm2.5 levels . existing methods used to estimate pm2.5 levels fall 
 in several classes that include ( 1 ) chemical transport models ( ctm ) , 
 ( 2 ) land use regression ( lur ) , ( 3 ) satellite data , and ( 4 ) different 
 geostatistical approaches . lur is a regression model which estimates 
 air pollution as a function of explanatory variables . lur takes characteristics 
 from the study area ( traffic count , road length , distance to nearest 
 road , elevation , land cover , household density , wind , etc . ) and develops 
 a multiple linear regression model which aims at describing a pollutant 
 of interest . most lur models are geared toward a model that explains 
 the most variability of the dependent variable ( i.e. , the model with 
 the highest possible r ) on a relatively 
 small scale . they range from process - based prediction methods to data - driven statistical 
 estimation methods . the first two classes of methods are defined by 
 their ability to predict levels based on a model representation of 
 the processes that lead to air pollution . the latter two classes 
 are driven by observations , such as satellite readings or measurement 
 from ground monitoring stations . although no 
 categorization is without exception or entirely distinctive , these 
 classes demonstrate possible methodological procedures . these four 
 classes also differ widely in terms of accuracy , complexity , numerical 
 cost and accessibility ( see supporting information 
 ( si ) ) . geostatistical approaches provide , like satellite 
 data , observationally 
 driven estimates of pm2.5 . they usually consist of linear 
 kriging estimators of pm2.5 calculated from measurements 
 at ground monitoring stations . these methods provide accurate estimates 
 in the vicinity of monitoring stations and are simple to use , thereby 
 providing a widely utilized approach . however , like any observationally 
 driven estimation method , geostatistical methods alone can not be used 
 to explain contributions from major contributing sources . while 
 work has been done to develop methods individually within 
 the four classes mentioned above , there is also interest in combining 
 approaches across classes to create an estimation framework that combines 
 the strengths of the respective groups . the goal of this work is to 
 combine a process based method and an observationally based estimation 
 method to create a combined estimation method that can be used by 
 a wide audience to accurately estimate the distribution the annual 
 pm2.5 concentration across the continental united states 
 ( u.s . ) from 1999 to 2009 , and to quantify how much of the estimated 
 annual pm2.5 concentration can be explained by the major 
 contributing sources of on road mobile emissions and stationary emissions . we will achieve our goal by using the bayesian maximum entropy 
 ( bme ) knowledge synthesis framework to combine 
 lur with geostatistical estimation . bme utilizes bayesian epistemic 
 knowledge blending to combine data from multiple sources . for our 
 process - based method we select lur over ctms because of its ability 
 to use readily available information about on road mobile emissions 
 and stationary emissions to predict annual pm2.5 . for our observationally based method we rely on a geostatistical analysis 
 of ground observations of pm2.5 concentrations because 
 of the relatively large number of monitoring stations providing accurate 
 measurements across the u.s . by combining methods like lur and bme 
 we can create a model that is numerically efficient , applicable and 
 interpretable over a large domain size . the knowledge base considered 
 in the bme method consists of general 
 knowledge describing generalizable characteristics of the space / time 
 pm2.5 field ( such as its space / time trends and dependencies , 
 its relationship with respect to various emissions , etc ) , and site 
 specific knowledge that include hard data ( data without measurement 
 error ) and soft data ( data with measurement errors which can be non - gaussian ) . the strategy we will use in this work is to employ lur to describe 
 the general trends of annual pm2.5 concentrations over 
 the entire u.s . and model the pm2.5 residuals ( obtained 
 by removing the lur offset ) using bme . this will allow us to rigorously 
 account for the non - gaussian uncertainty associated with annual pm2.5 concentration calculated from daily concentrations where 
 some daily concentrations may be missing . one outcome of our work is the development of an lur for the 
 prediction 
 of annual pm2.5 concentrations across the continental u.s . , 
 which is a geographical domain of a fairly large size . while many 
 previous studies have developed lur models over small geographical 
 domains where high predictability can be achieved , each specific lur model is usually only valid for the small 
 region for which it was developed . in 
 other words there have been comparatively 
 fewer studies that developed an lur with lower predictability but 
 higher generalizability . the lur we present fills that knowledge gap , 
 with a specific focus on using annual pm2.5 explained by 
 on road mobile emissions and stationary emissions as its predictors . another outcome of our work is the sequential integration of two 
 classes of methods ( lur and geostatistical ) to create a combined lur / bme 
 estimation method that borrows strengths from each of its constituent . combining methods is a growing research area and our work contributes 
 to that field . while very few works have looked at combining lur and 
 bme approaches or lur and kriging approaches , more studies are needed in order to explore 
 the various ways by which to combine these methods . we focus specifically 
 on using lur to provide general knowledge about pm2.5 , 
 using bme to account for the incompleteness of daily samples , and 
 making the combined method accessible to a wide audience . other strategies 
 and focus will undoubtedly have to be investigated in future works , 
 for example creating more elaborate lur models including those which use meteorological data . finally we use our lur / bme model to perform 
 a risk assessment that 
 differentiates the number of annual pm2.5 predicted deaths 
 that can be explained by on road mobile emissions and stationary emissions . the dichotomous assignment of pm2.5 to these two sources 
 allows for straightforward abatement strategies . this assessment is 
 useful on its own to generate research questions that can improve 
 methods used to calculate death reductions achieved under various 
 scenarios of source reductions . raw daily federally referenced 
 method ( frm ) pm2.5 monitoring data collected from 1999 
 to 2009 were obtained from the epa s air quality systems ( aqs ) 
 database across the contiguous united states . whenever a daily pm2.5 monitoring value reported below 
 the detection limit of its monitor , it was replaced by the mean of 
 a log - normal distribution that was fit to all reported below - detect 
 values . daily values were averaged whenever two or more daily pm2.5 monitoring values were reported by collocated monitors 
 on a given day / site . annual pm2.5 were calculated 
 from daily pm2.5 monitoring values as follows : every day 
 for which a station reported a daily pm2.5 monitoring value , 
 a corresponding annual pm2.5 was calculated by taking the 
 arithmetic average of all the daily monitoring values reported at 
 that station over the previous year ( i.e. , 365 days ) including that 
 day . note that this one year period could include time before january 
 1 , 1999 ( i.e. , the first day for which daily monitoring data were 
 available ) . the intended sampling frequency of a given daily 
 monitoring station 
 was used to calculate how many daily monitoring values should have 
 been reported in a given year period . comparing this number to the 
 actual number of reported monitoring values informs us about the incompleteness 
 of intended sampling over that given year . the lur model predicts annual pm2.5 given 
 a series of predictive lur independent variables that characterize 
 the effect of ( a ) elevation , ( b ) on road mobile emissions , and ( c ) 
 stationary emissions . a detailed explanation of all data sources for 
 each lur independent variable is described in the si . we focus on on road mobile emissions and stationary 
 emissions because they are two major contributors to anthropogenic 
 pollution . for stationary emissions , we used data from the epa s 
 national emissions inventory ( nei ) , which 
 provides inventories of stationary emissions ( in tons / year ) of the 
 main constituents of pm2.5 ( i.e. , so2 , nh3 , pm2.5-primary and nox ) . we assume that at space / time location p = ( s , t ) , the effects of stationary source emissions decrease exponentially 
 with distance between the source and the location p , as given by the equation vi , p = n  =  1nemi(en , t)exp((3en p)/(dri ) ) , where , i = so2 , nh3 , pm2.5 , nox , emi(en , t ) is the emissions in tons / year of constituent i at stationary source emissions location en and time t and dri is the exponential decay range in km . it would be difficult to consistently 
 and accurately measure on road mobile emissions across the entire 
 u.s . thus for on road mobile emissions we use data estimating vehicular 
 traffic ( annual average daily traffic counts for each major highway 
 road segment in the u.s . as estimated through linear referencing ) and population density ( people / km ) to construct variables that estimate total traffic ( tt ) , average 
 congestion ( ac ) , and emission efficiency ( ee ) based on population 
 density . emission efficiency is added to correct for the assumption 
 that every mile driven produces the same amount of emissions regardless 
 of vehicle type by hypothesizing that areas with high population density 
 tend to have vehicles better suited for urbanized environments , which 
 ( in general ) are more fuel efficient . these traffic and emission efficiency 
 variables are then combined to provide an estimate of on road mobile 
 emission , thereby bypassing the laborious task of obtaining on road 
 mobile emission data directly for a nationally sized domain . our large area variability 
 lur expresses the annual pm2.5 at space / time location p = ( s , t ) , where s = ( s1,s2 ) is the spatial coordinate 
 and t is time , as a linear combination of the corresponding 
 lur independent variables at p. the 
 first independent variable consists of the elevation velevv , p at p. the next three independent 
 variables characterize the effect of on road mobile emissions . they 
 are denoted as the column vector vmobile , p = [ vtt , pvac , pvee , p ] , where the subscript t denotes the transpose , and vtt , p , vac , p and vee , p are variables characterizing total traffic , average 
 congestion , and emission efficiency , respectively , at p. the last four independent variables characterize 
 the effect of stationary emissions . they are denoted as vstationary , p = [ vso2,pvnh3,pvpm2.5,pvnox , p ] , where vso2,p , vnh3,p , vpm2.5,p and vnox , p are variables charactering the concentrations of so2 , nh3 , pm2.5 , and nox , respectively at space / time location p. we consider models that include the elevation variable , 
 at least 1 out of the 3 on road mobile emission variables , and at 
 least 1 out of the 4 stationary emission variables , which results 
 in a total of   candidate models . these models are expressed 
 by the following equation1where zp is annual pm2.5 at p , 0 is the equation 
 intercept , 
 elev , mobile = [ tt ac ee ] and stationary = [ so2 nh3 pm2.5 nox ] are linear coefficients 
 for the independent variables velev , p , vmobile , p and vstationary , p , respectively , imobile = [ ittiaciee ] and istationary = [ iso2inh3ipm2.5inox ] are vectors of indicator values ( 0 or 1 ) such that at 
 least one element in both imobile and ipoint must be 1 , the operator 
 denotes the element - by - element multiplication between same - sized vectors 
 and p is a homoscedastic 
 error term . due to the large overlap in annual pm2.5 , only a subset 
 of annual pm2.5 was used to construct the lur model to 
 avoid collinearity . namely , only the last annual pm2.5 in 
 a calendar year was used from each station ( approximately 11  000 
 data values ) , encompassing all daily values . each of the 105 
 candidate lur models were optimized by selecting 
 hyperparameter values that maximized the lur r. a hyperparameter is a physical parameter within each variable 
 that is allowed to adjust based on predictability of annual pm2.5 . hyperparameters for annual pm2.5 include the 
 radii ar1 , ar2 and ar3 for the buffers used to 
 calculate total traffic , average congestion , and emission efficiency , 
 respectively , and the exponential decay ranges for stationary source 
 variables ( i.e. , dri described 
 in the si ) . the fminsearch function of matlab was used to search for hyperparameter values 
 that maximized the lur r. the search 
 was started given an initial selection of hyperparameters described 
 in si . the akaike information criteria 
 ( aic ) and all variance inflation 
 factor ( vif ) values were found for each of the 105 optimized candidate 
 lur models . aic is a measure of parsimony of a model and vif is a 
 measure of collinearity of a model . out of the 105 optimized models , 
 our final model has the lowest aic value among models with vif values 
 < 10 and with physically plausible s . the s have to 
 be positive in order to be plausible , with the exception of negative 
 s for emission efficiency and elevation . bme can incorporate information from many different sources and 
 bme is implemented using the bmelib suite of functions 
 in matlab . the buttress of bme has been 
 detailed in other works , and can be summarized as performing the following steps : ( 1 ) gathering 
 the general knowledge base ( g - kb ) and site - specific knowledge base 
 ( s - kb ) about the mapping situation , ( 2 ) using the maximum 
 entropy principle of information theory to process the g - kb 
 in the form of a prior probability distribution function ( pdf ) fg , ( 3 ) integrating s - kb using 
 an epistemic bayesian conditionalization rule on 
 data fs with and without 
 measurement error to create a posterior pdf fk , and ( 4 ) creating space / time estimates based 
 on the analysis . we use a space / time random field ( s / trf ) to describe 
 the variability of annual pm2.5 across the u.s . our notation 
 a for s / trf will consist of denoting a single random variable z in capital letters , its realization , z , in lower case ; and vectors and matrices in bold faces ( e.g. , z = [ z1, ... ,zn ] and z = [ z1, ... ,zn ] ) . let z(p ) = z(s , t ) be a space / time random field ( s / trf ) representing annual pm2.5 . we define the transformation of the pm2.5 data zh observed at locations ph as2where oz(p ) may be any deterministic 
 offset that can be calculated without error as a function of the space / time 
 coordinate p. we then define x(p ) as the s / trf representing 
 the variability and uncertainty associated with the transformed data xh , and we let z(p ) = x(p)+oz(p ) be the s / trf representing 
 pm2.5 . in this work , we consider two choices for oz(p ) : ( 1 ) a 
 constant value and ( 2 ) the lur estimate lur , p given by3where the estimated s indicators and s coefficients are those derived in 
 our final annual pm2.5 lur model . we can then calculate k , the estimated annual 
 pm2.5 at unmonitored location pk by obtaining the bme estimate xk for the transformed 
 s / trf x(p ) at the estimation 
 point pk , and adding back oz(pk ) , the 
 offset calculated at pk . the g - kb for the transformed s / trf x(p ) consists of its expected value mx(p ) and 
 covariance function cx(p , p ) 
 ( see si ) . the s - kb for x(p ) consists of hard and soft data . 
 the hard data xh = zh  oz(ph ) are obtained 
 based on annual pm2.5 values zh calculated at 
 hard data points ph where at least 75% of intended samples were collected , in 
 line with epa regulations pertaining to valid design values . data points not meeting this completeness criterion 
 are classified as the soft data points ps , with an uncertainty attributed to 
 the incompleteness of intended sampling . following akita et al , the uncertainty associated with the annual pm2.5 , zs for station i and date t is described by a gaussian 
 pdf truncated below zero , with mean i , t and standard deviation i , t . the mean i , t is simply the sample mean of 
 the ni , t daily concentrations ( zi , j , j = 1, ... ,ni , t ) recorded at station i over 1 year preceding date t. the epistemic 
 uncertainty associated with the incompleteness of intended sampling 
 is characterized by the difference between ni , t and the intended number 
 of samples ni *  ni , t that would 
 have been collected if the station worked as intended in accordance 
 with the monitor s sampling frequency . therefore a reasonable 
 choice for the standard deviation quantifying that uncertainty is4where the first factor is the standard deviation 
 of the sample mean and the second factor is a population correction 
 factor that accounts for the incompleteness of intended sampling from 
 a population of size ni*. the pdf for xs is 
 then derived from the pdf for zs by simply using the transformation xs = zs the 
 g - kb and s - kb for the s / trf x(p ) can overall be written as g = 
 { mx(p),cx(p , p ) } 
 and s = { xh , fs ( xs ) } , 
 and in this case the bme posterior pdf for x(pk ) at estimation 
 point pk is given by fk(xk ) = a  dxfs(xs ) fg(x ) where x = ( xk , xh , xs ) is a realization of x at points p = ( pk , ph , ps ) and a is a normalization constant . finally the pdf for zk is obtained by simply using the back - transformation zk = xk + oz(pk ) . in order to test 
 the estimation improvement of lur and bme , a cross - validation 
 was performed to compare three different methods used in this study : 
 ( a - constant / hard ) setting the deterministic global 
 ( i.e. , covering a substantial domain where variability within the 
 domain can be largely diverse ) offset oz(p ) to a constant value 
 and considering all data as hard , ( b - lur / hard ) setting 
 the global offset to the lur model and considering all data as hard 
 and ( c - lur / hard and soft ) setting the global offset 
 to the lur predicted value and considering data as hard and soft as 
 defined in the previous section . for each of these methods , the cross 
 validation procedure consists of randomly selecting 20  000 
 hard data points , removing each one at a time , and re - estimating it 
 from the remaining annual pm2.5 . the cross - validation statistics 
 investigated include mean squared estimation error ( mse ) , root mean 
 squared estimation error ( rmse ) , mean absolute estimation error ( mae ) , 
 mean of the root variance of the posterior pdfs ( mr ) , the square of 
 pearson s correlation coefficient , and the square of spearman s 
 correlation coefficient . along with the leave - one - out cross 
 validation ( loocv ) of 20 ,  000 hard data point , a 10-fold spatial 
 cross - validation was also performed . the incorporation of the 
 lur model into the bme methodology has many potential applications 
 including determining the mortality of various diseases attributable 
 to pm2.5 . excess mortality was calculated using the methodology 
 presented by li , assuming linearity , 
 in order to quantify total mortality , mortality from ischemic heart 
 disease ( ihd ) and mortality from lung cancer ( lc ) . excess deaths were calculated for ( 1 ) annual 
 pm2.5 , ( 2 ) annual pm2.5 explained by on road 
 mobile emissions , and ( 3 ) annual pm2.5 explained by stationary 
 emissions . let l ( p ) denote our estimate of annual 
 concentrations , where l = total for total pm2.5 , l = mobile for pm2.5 explained 
 by on road mobile emissions , and l = stationary for 
 pm2.5 explained by stationary emissions . for l = total we simply use total ( p ) = lur / bme , p , where lur / bme,(p ) is the lur / bme 
 estimate of annual pm2.5 described earlier . for l = mobile we use the lur in a relative manner to estimate 
 the ratio lur  mobile(p ) = ( ( mobile  mobile)vmobile , p)/(lur , p ) corresponding to the proportion of pm2.5 that 
 the lur model explains from on road mobile emissions . we then multiply 
 that ratio with the lur / bme estimate of annual pm2.5 , so 
 that mobile(p ) = lur / bme , plur  mobile(p ) . likewise we use stationary(p ) = lur / bme , plur  stationary(p , with lur  stationary(p ) = ( ( stationary  stationary)vstationary , p)/(lur , p ) ) . the mortality for a specific cause 
 of death ( e.g. , total mortality , ihd , lc ) attributed to an annual 
 concentration l ( p ) is given by li et al.5where i0 is the 
 baseline incidence rate for the cause of death of interest ,  
 is the corresponding concentration response coefficient , p is the population at the county level , and zb is the background concentration . for pm2.5 of 35 g / m . we use zb = 5 g / m . there were 1  478  149 
 annual pm2.5 data points from 1999 to 2009 coming from 
 1576 monitoring stations . the mean of the annual pm2.5 is 12.44 g / m , the variance is 11.57(g / m ) , the skewness is 0.56 and the kurtosis is 5.57 . the minimum annual 
 value is 1.63 g / m and the maximum annual value 
 is 75.40 g / m . the final lur model 
 had six independent variables : elevation , three on road mobile emission 
 variables ( total traffic , average congestion , emission efficiency ) , 
 and two stationary emission variables ( nh3 and so2 ) ( table 1 ) . table 1 describes the optimal hyperparameters for each variable along with 
 their corresponding  values . this lur model has an r = 0.53 , providing generalizable predictability 
 of annual pm2.5 over the entire u.s . from 1999 to 2009 . the combination of 
 the lur and bme methods 
 through methods ( a ) to ( c ) led to a refined estimation of annual pm2.5 as seen in figure 1 showing estimated 
 levels across the u.s . for may 1 , 1999 . method ( a - constant / hard ) using a constant offset and using all data as hard does not differentiate 
 well the annual pm2.5 across southern california and estimates 
 fairly benign levels for several states west of the mississippi river . 
 by incorporating the lur offset , method ( b - lur / hard ) provides estimates of annual pm2.5 that are more refined 
 and localized . by further incorporating the soft data to the hard 
 data and lur offset , method ( c - lur / hard and soft ) 
 further refines the description of hot spots across the country . method 
 ( c ) is able to pick up finer scale variation in concentrations compared 
 to methods ( a ) and ( b ) . this finer scale variation can also be seen 
 in subsequent months ( si figure s5 ) . bme predicted 
 annual pm2.5 ( g / m ) 
 concentration estimation map across the contiguous u.s . on may 1 , 
 1999 for the following methods : ( a ) constant offset / hard data ; ( b ) 
 lur offset / hard data ; and ( c ) lur offset / hard and soft data . cross validation statistical measures 
 indicated a consistent improvement 
 in mapping accuracy from method ( a ) to ( c ) ( table 2 ) . measures of estimations errors ( mse , rmse , mae , mr ) decreased 
 from method ( a ) to ( b ) and from method ( b ) to ( c ) , while measures 
 of correlation ( square pearson s corr . . coeff . ) increased from method ( a ) to ( b ) and from method ( b ) 
 to ( c ) . incorporating the lur offset while using only hard data ( i.e. , 
 going from method ( a ) to ( b ) ) resulted in a reduction of 21.89% in 
 mse . further incorporating soft data ( i.e. , going from method ( b ) 
 to ( c ) ) resulted in an additional reduction of 4.87% in mse . the reduction 
 in mse from method ( b ) to ( c ) is more pronounced when performing cross - validation 
 on points that contain a higher percentage of soft data ( si table s3 ) . this reduction is more pronounced 
 still when estimation neighborhoods around cross - validation locations 
 are forced to have soft data points ( si table s4 ) . this corresponds to 12.8% shrinkage in r , which is reasonable since the training set for the 
 10-fold cross validation is substantially smaller than that of the 
 loocv . using eq 5 with total(p ) we find that the number of deaths from 1999 to 2007 predicted from 
 annual pm2.5 exposure in excess of background levels is 
 905  560 . these results were validated using the epa s 
 benmap program and are consistent with 
 other estimates . we then used eq 5 with mobile(p ) ( pm2.5 explained by on road mobile 
 emissions ) and stationary(p ) ( pm2.5 explained by stationary 
 emissions ) . the mean of the mobile(p ) across all the space / time data 
 points is 3.4 g / m , while the mean of stationary(p ) across the 
 same points is only 1.15 g / m . accordingly the number 
 of deaths attributed to pm2.5 explained by on road mobile 
 emissions is greater than the number of deaths attributed to pm2.5 explained by stationary emissions ( table 3 ) . for instance , the number of deaths attributed to pm2.5 
 explained by on road mobile emissions is 568  090 from 1999 
 to 2007 , which is 1.85 times more than the 306  316 deaths attributed 
 to pm2.5 explained by stationary emissions . similarly , 
 on road mobile emissions explained 1.86 times the number of ihd deaths 
 and 1.98 times the number of lc deaths compared to deaths explained 
 by stationary emissions . the number of deaths assumes that the relative 
 risk used in eq 5 can be applied to the entire 
 population and that estimated ambient concentration is a surrogate 
 for exposure . this finding is interesting 
 because , according to the nei , primary 
 pm2.5 , no2 , so2 , and nh3 coming from on road mobile emissions sum up to 70  834 thousand 
 tons from 1999 to 2007 while primary pm2.5 , no2 , so2 , and nh3 coming from stationary emissions 
 sum up to 293  446 thousand tons for the same time period ( si table s2 ) . hence , even though on road mobile 
 emissions emit only about a quarter of the mass emitted by stationary 
 emissions , the number of deaths predicted from pm2.5 explained 
 by on road mobile emissions is almost twice that predicted from pm2.5 explained by stationary emissions . the first major outcome of our work is the creation of a global 
 lur model that predicts large area variability of pm2.5 across the entire contiguous united states from 1999 to 2009 . only 
 a handful of studies have developed lur models that can be classified 
 as  general  in that they produced results generalizable 
 to domain sizes as large as ours ( si figure 
 s1 ) . although the lur may perform better in some areas than others , 
 the model is generalizable  in a relative fashion when 
 compare to lur models developed over a smaller domain . to the best 
 of our knowledge , the closest lur models developed over such a large 
 domain size are hart et al . and beelen 
 et al . for annual pm10 , and 
 beckerman et al . for monthly pm2.5 . studies developed regression models 
 to predict 
 annual pm10 concentrations across the united states from 
 1985 to 2000 and across 15 european countries for 2001 , respectively . model used traffic related variables while the 
 beelen et al . model used meteorology 
 and land use ) , they produced similar r of 49% and 41% , respectively . our study is successful in helping to fill that knowledge gap by providing 
 a general lur for pm2.5 that achieves an r of 53% that is comparable or better than that for annual 
 pm10 . of the limited general lur models developed 
 for the long - term average 
 concentration of pm2.5 , the lur - without - remote - sensing 
 model developed by beckerman et al . the explanatory variables of that model 
 are traffic within 1 km and green space within 0.1 km . the r of that model was 3% for their training data set and 
 5% for their validation data set . this provides a substantial contribution 
 to the literature because it describes the small area variability 
 at a scale of 0.11 km which is attributable to primary pm2.5 emitted as ultrafine particles and specific metals that 
 are responsible for a large portion of observed health impacts . hence their relatively small r of 35% makes sense since on road emission of 
 primary pm2.5 makes up only a small fraction of annual 
 pm2.5 . by contrast , our explanatory variables capture the 
 large area variability of pm2.5 , as demonstrated by our 
 hyperparameter values ranging from ten to hundreds of km ( table 1 ) . these values are characteristic of the transport 
 of secondary pm2.5 over long distances . since the majority of pm2.5 comes from secondary formation , 
 it makes sense that our model achieves a higher r ( 53% ) than that of beckerman et al . therefore their model describes short - area 
 variability due to the local effect of traffic coming from primary 
 pm2.5 , while our model describes the large area variability 
 of pm2.5 due to the long - range effect of secondary pm2.5 created from on road mobile emissions and stationary emissions . 
 to the best of our knowledge no other lur models have predicted the 
 large area variability of secondary pm2.5 and our model 
 is the first to address this important knowledge gap . our lur model 
 estimated pm2.5 coming from on road mobile emissions using 
 tt , ac , and ee . however , due to the complex nature of pm2.5 , these variables may capture other sources outside of on road mobile 
 emissions . we hypothesized that ee using population density corrects 
 for over prediction of on road mobile emissions coming from tt and 
 ac . indeed , out of the   models that had the ee variable , ee was positive for the 1  ( 1 + 1 + 
 0 )  15 = 30 models where ee appears without tt and it consistently 
 switched to being negative for the 1  ( 0 + 1 + 1 )  15 
 = 30 models that contain both the ee and the tt variable . however , 
 when paired with the tt traffic variable , ee corrects the overestimation 
 of these variables . this finding is in agreement of our hypothesis 
 and therefore supports the conclusion that population density can 
 be used as a surrogate for increased ee of the vehicle fleet . obtaining 
 accurate estimates of on road mobile emissions along all roads is 
 a difficult task . by using population data to calculate ee , we facilitate 
 this task and as a result we ensure the accessibility of our model 
 to a wider audience . while previous lur models represent important 
 contributions to 
 the field , our model differs in several important ways : ( 1 ) our model 
 describes large area variability of pm2.5 , which characterizes 
 the secondary component of this pollutant , ( 2 ) the explanatory variables 
 are constructed from data that are easily obtainable by a wide audience 
 and ( 3 ) our model allows to distinguish between pm2.5 explained 
 by on road mobile emissions and pm2.5 explained by stationary 
 emissions . to our knowledge this is one of the first lur models to 
 capture secondary pm2.5 using easily obtainable explanatory 
 variables describing on road mobile emissions and stationary emissions . the second major outcome of this work is the combination of our 
 lur model with bme to create a combined lur / bme hybrid estimation 
 method for annual pm2.5 . in this hybrid approach , lur is 
 used as a first step to characterize global trends in pm2.5 and bme is used to extract unexplained variability in the residuals . 
 our results ( table 2 ) demonstrate that lur / bme 
 is successful at combining the strengths of each of its component 
 methods . indeed , lur / bme results in a 21.89% reduction in mse and 
 a 28.94% increase in r over bme alone , 
 which is itself more accurate than lur alone . the population correction 
 factor presented in the soft data variance in eq 4 does not account for the fact that annual pm2.5 averages 
 are correlated in time . as well , the number of daily values within 
 a year ni , t does not account for the seasonality of missing values . their work 
 saw an r of 0.79 using a validation data 
 set comprised of about 10% of the data . by comparison we achieved 
 an r of 0.78 using a 10-fold cross validation , 
 where each of validation points had similar distance - to - closest - monitor 
 as those of beckerman et al . a key difference between our works is 
 that we extended their work by incorporating non - gaussian soft data 
 that rigorously accounted for the uncertainty associated with the 
 incompleteness of daily samples . our r indicates that our model was successful in this novel incorporation 
 of non - gaussian soft data in the lur / bme framework , which resulted 
 in one of the most accurate lur / bme estimations to date of annual 
 pm2.5 as supported by the fact that our r is similar to that of beckerman et al . a unique strength 
 of our model is that these highly accurate lur / bme estimates of annual 
 pm2.5 can be separated into the portions explained by on 
 road mobile emissions and stationary emissions , which to our knowledge 
 had not been done before to a similar level of precision . building 
 on the novel contributions of the first two outcomes of 
 our work , an important third outcome of this work is a risk assessment 
 of annual pm2.5 exposure explained from major contributing 
 sources . estimating annual pm2.5 is useful for assessing 
 long - term exposure needed to investigate chronic diseases . others 
 have already used lur estimates in epidemiological studies . from 1999 to 2007 there were 568  090 deaths attributed to 
 pm2.5 explained by 70  834 thousand tons of primary 
 pm2.5 , no2 , so2 , and nh3 emitted 
 by on road mobile emissions , which correspond to a ratio of 8.02 deaths / thousand 
 tons for on road mobile emissions . by contrast there were 306  316 
 deaths attributed pm2.5 explained by 293  446 thousand 
 tons of primary pm2.5 , no2 , so2 , 
 and nh3 emitted by stationary emissions , which correspond 
 to a ratio of 1.04 deaths / thousand tons for stationary emissions . these results are informative because they imply that mechanisms involved 
 in the creation and long - range transport of secondary pm2.5 lead to substantially differing health impacts depending on whether 
 emissions originate from on road mobile emissions or stationary emissions . other works have also examined excess mortality due to current 
 emissions levels . when investigating massachusetts power plants levy 
 and spengler found that current power 
 plant emissions in the surrounding area that emitted above the best 
 available control technology ( bact ) resulted in approximately 70 deaths 
 per year in a  600 km by 600 km region which includes areas 
 of massachusetts and new york where the power plants were located . according to the bact of 3 lb / mwh of so2 and 1.5 lb / mwh 
 of nox , there would be a reduction of 
 43  951 tons of so2 and 4376 tons of nox from the two power plants mentioned in the study . this would result in 1.34 deaths / thousand tons of so2 and 2.51 deaths / thousand tons of nox due to power plants emissions in the area being above the 
 bact . that work used the ctm calpuff in which emission levels can 
 be zeroed out while our work uses an lur model which measures annual 
 predicted pm2.5 . even though lur can not be directly 
 compared to ctms , our lur results are useful in a relative manner 
 as they allow us to contrast on road mobile emissions and stationary 
 emissions which have not been done before . in order to reduce 
 the number of deaths due to pm2.5 exposure , our results 
 indicate a reduction in one ton of on road 
 mobile emissions would be eight times more beneficial than a one ton 
 reduction in stationary emissions . this may be accomplished though 
 any number of actions such as increased accessibility and reliance 
 on public transportation in areas of high population density to more 
 stringent emission standards that would further promote fuel efficiency .
knowledge 
 of particulate matter concentrations < 2.5 m 
 in diameter ( pm2.5 ) across the united states is limited 
 due to sparse monitoring across space and time . 
 epidemiological studies 
 need accurate exposure estimates in order to properly investigate 
 potential morbidity and mortality . 
 previous works have used geostatistics 
 and land use regression ( lur ) separately to quantify exposure . 
 this 
 work combines both methods by incorporating a large area variability 
 lur model that accounts for on road mobile emissions and stationary 
 source emissions along with data that take into account incompleteness 
 of pm2.5 monitors into the modern geostatistical bayesian 
 maximum entropy ( bme ) framework to estimate pm2.5 across 
 the united states from 1999 to 2009 . 
 a cross - validation was done to 
 determine the improvement of the estimate due to the lur incorporation 
 into bme . 
 these results were applied to known diseases to determine 
 predicted mortality coming from total pm2.5 as well as 
 pm2.5 explained by major contributing sources . 
 this method 
 showed a mean squared error reduction of over 21.89% oversimple kriging . 
 pm2.5 explained by on road mobile emissions and stationary 
 emissions contributed to nearly 568  090 and 306  316 
 deaths , respectively , across the united states from 1999 to 2007 .