Code for Bayesian Hierarchical Spatial Modeling of New York City ZIP Code-Level Positive COVID19 Testing Results
Charles DiMaggio, PhD, MPH
Bellevue-NYU Department of Surgery
Division of Trauma and Critical Care
Last update: 04 May, 2020
Abstract
Introduction: The population and spatial characteristics of COVID-19 infections are poorly understood.
Methods: We analyzed positive COVID19 testing results counts within ZIP Code Tabulation Areas (ZCTA) with Bayesian hierarchical Poisson spatial models using integrated nested Laplace approximations.
Results: Spatial clustering accounted for approximately 32% of the variation in the data, with hot spots in all five boroughs. Spatial risk did not correspond precisely to population-based rates of positive tests. The strongest univariate association was with the proportion of residents in a ZIP Code Tabulation Area with Chronic Obstructive Pulmonary Disese (COPD). For every one unit increase in a scaled standardized measure of COPD, there was an approximate 8-fold increase in the risk of a positive COVID-19 test in a ZCTA (Incidence Density Ratio = 8.2, 95% Credible Interval 3.7, 18.3). The next strongest association with with the proportion of Black residents, for which there there was a nearly five-fold increase in the risk of a positive COVID19 test. (IDR = 4.8, 95% Cr I 2.4, 9.7). The proportion of residents older than 65, housing density and the proportion of residents with heart disease were each associated with an approximate doubling of risk. In a multivariable model including estimates for age, COPD, heart disease, housing density and black race, the only variables that remained associated with positive COVID-19 testing with a probability greater than chance were the proportion of black residents and older persons.
Conclusions: The population and spatial patterns of COVID-19 infections differ by race, age, physical environment and health status. Areas with large proportions of Black residents are at markedly higher risk that is not explained by characteristics of the environment and pre-existing conditions in the population.
1 Load Data and Packages
Testing data were obtained from the New York City Department of Health and Mental Hygeine (NYC DOHMH) releasing COVID19 tests results posted on GitHub. A detailed description of data creation can be found below. Covariate are drawn primarily from 2010 US Census results.
Prep the data by ordering to match the spatial dataframe, and add a sequential area ID variable for the model statements.
library(maptools, quietly = T)
library(sp)
library(spdep, quietly = T)
library(INLA, quietly = T)
library(ggplot2)
library(kableExtra)
library(DescTools)
## Registered S3 method overwritten by 'DescTools':
## method from
## reorder.factor gdata
options(scipen = 10) # don't use scientific notation
options(fmt.num=structure(list(digits=2, big.mark=","), class="fmt")) # force use commas as thousands separator and 2 decimal places
testDat<-readRDS("~/testDat.RDS")
nyc <- readRDS("~/nycMap4COVID.RDS")
nycDat <- attr(nyc, "data")
order <- match(nycDat$ZCTA5CE10,testDat$zcta)
testDat<- testDat[order,]
testDat$ID.area<-seq(1,nrow(testDat))
nyc.adj <- paste("~/nyc.graph")
1.1 Descriptive Statistics
There were 177 ZCTA’s in the data set. The mean COVID-19 test rate per 10,000 ZCTA population was 166.2 (95% CI 156.7, 175.7). The distribution appeared multimodal. The mean COVID-19 test rate per 10,000 tests was 5,176.0 (95% CI 5,045.9, 5,306.1) and appeared skewed. The The 5 ZCTAs with highest positive COVID-19 test numbers per 10,000 population were the same as those with the highest proportion per 10,000 tests (10464, 10470, 10455, 10473, 11234, and 11210). The 5 lowest ZCTAs were also the same for both measures 1(1103, 11102, 11693, 11369, 11363, and 10308). Table one presents comparative statistics for the ZCTA’s with the highest and lowest quantiles for population-based positive test rates.
## -------------------------------------------------------------------------
## Positive COVID-19 Tests per 10,000 ZCTA Population
##
## length n NAs unique 0s mean
## 177 177 0 = n 0 166.18312
## 100.0% 0.0% 0.0%
##
## .05 .10 .25 median .75 .90
## 71.10536 84.17626 113.19039 162.97451 215.08084 248.76561
##
## range sd vcoef mad IQR skew
## 281.17539 64.25910 0.38668 75.00825 101.89045 0.18618
##
## meanCI
## 156.65093
## 175.71531
##
## .95
## 267.41561
##
## kurt
## -0.81984
##
## lowest : 42.02216, 49.81734, 50.22831, 60.10303, 61.66783
## highest: 282.27373, 297.51712, 316.63043, 322.37247, 323.19755
## -------------------------------------------------------------------------
## Positive COVID-19 Tests per 10,000 ZCTA Tests
##
## length n NAs unique 0s mean
## 177 177 0 = n 0 5'175.967
## 100.0% 0.0% 0.0%
##
## .05 .10 .25 median .75 .90
## 3'608.459 3'907.126 4'519.016 5'380.117 5'852.102 6'039.269
##
## range sd vcoef mad IQR skew
## 4'747.624 876.960 0.169 780.456 1'333.086 -0.589
##
## meanCI
## 5'045.879
## 5'306.055
##
## .95
## 6'198.431
##
## kurt
## -0.128
##
## lowest : 2'586.207, 2'937.063, 2'972.973, 3'142.857, 3'178.295
## highest: 6'481.481, 6'614.493, 6'818.530, 6'982.703, 7'333.831
Charles DiMaggio, PhD, MPH
Bellevue-NYU Department of SurgeryDivision of Trauma and Critical Care
Last update: 04 May, 2020
Abstract
Introduction: The population and spatial characteristics of COVID-19 infections are poorly understood.Methods: We analyzed positive COVID19 testing results counts within ZIP Code Tabulation Areas (ZCTA) with Bayesian hierarchical Poisson spatial models using integrated nested Laplace approximations.
Results: Spatial clustering accounted for approximately 32% of the variation in the data, with hot spots in all five boroughs. Spatial risk did not correspond precisely to population-based rates of positive tests. The strongest univariate association was with the proportion of residents in a ZIP Code Tabulation Area with Chronic Obstructive Pulmonary Disese (COPD). For every one unit increase in a scaled standardized measure of COPD, there was an approximate 8-fold increase in the risk of a positive COVID-19 test in a ZCTA (Incidence Density Ratio = 8.2, 95% Credible Interval 3.7, 18.3). The next strongest association with with the proportion of Black residents, for which there there was a nearly five-fold increase in the risk of a positive COVID19 test. (IDR = 4.8, 95% Cr I 2.4, 9.7). The proportion of residents older than 65, housing density and the proportion of residents with heart disease were each associated with an approximate doubling of risk. In a multivariable model including estimates for age, COPD, heart disease, housing density and black race, the only variables that remained associated with positive COVID-19 testing with a probability greater than chance were the proportion of black residents and older persons.
Conclusions: The population and spatial patterns of COVID-19 infections differ by race, age, physical environment and health status. Areas with large proportions of Black residents are at markedly higher risk that is not explained by characteristics of the environment and pre-existing conditions in the population.
1 Load Data and Packages
Testing data were obtained from the New York City Department of Health and Mental Hygeine (NYC DOHMH) releasing COVID19 tests results posted on GitHub. A detailed description of data creation can be found below. Covariate are drawn primarily from 2010 US Census results.
Prep the data by ordering to match the spatial dataframe, and add a sequential area ID variable for the model statements.
library(maptools, quietly = T)
library(sp)
library(spdep, quietly = T)
library(INLA, quietly = T)
library(ggplot2)
library(kableExtra)
library(DescTools)## Registered S3 method overwritten by 'DescTools':
## method from
## reorder.factor gdataoptions(scipen = 10) # don't use scientific notation
options(fmt.num=structure(list(digits=2, big.mark=","), class="fmt")) # force use commas as thousands separator and 2 decimal places
testDat<-readRDS("~/testDat.RDS")
nyc <- readRDS("~/nycMap4COVID.RDS")
nycDat <- attr(nyc, "data")
order <- match(nycDat$ZCTA5CE10,testDat$zcta)
testDat<- testDat[order,]
testDat$ID.area<-seq(1,nrow(testDat))
nyc.adj <- paste("~/nyc.graph")1.1 Descriptive Statistics
There were 177 ZCTA’s in the data set. The mean COVID-19 test rate per 10,000 ZCTA population was 166.2 (95% CI 156.7, 175.7). The distribution appeared multimodal. The mean COVID-19 test rate per 10,000 tests was 5,176.0 (95% CI 5,045.9, 5,306.1) and appeared skewed. The The 5 ZCTAs with highest positive COVID-19 test numbers per 10,000 population were the same as those with the highest proportion per 10,000 tests (10464, 10470, 10455, 10473, 11234, and 11210). The 5 lowest ZCTAs were also the same for both measures 1(1103, 11102, 11693, 11369, 11363, and 10308). Table one presents comparative statistics for the ZCTA’s with the highest and lowest quantiles for population-based positive test rates.
## -------------------------------------------------------------------------
## Positive COVID-19 Tests per 10,000 ZCTA Population
##
## length n NAs unique 0s mean
## 177 177 0 = n 0 166.18312
## 100.0% 0.0% 0.0%
##
## .05 .10 .25 median .75 .90
## 71.10536 84.17626 113.19039 162.97451 215.08084 248.76561
##
## range sd vcoef mad IQR skew
## 281.17539 64.25910 0.38668 75.00825 101.89045 0.18618
##
## meanCI
## 156.65093
## 175.71531
##
## .95
## 267.41561
##
## kurt
## -0.81984
##
## lowest : 42.02216, 49.81734, 50.22831, 60.10303, 61.66783
## highest: 282.27373, 297.51712, 316.63043, 322.37247, 323.19755
## -------------------------------------------------------------------------
## Positive COVID-19 Tests per 10,000 ZCTA Tests
##
## length n NAs unique 0s mean
## 177 177 0 = n 0 5'175.967
## 100.0% 0.0% 0.0%
##
## .05 .10 .25 median .75 .90
## 3'608.459 3'907.126 4'519.016 5'380.117 5'852.102 6'039.269
##
## range sd vcoef mad IQR skew
## 4'747.624 876.960 0.169 780.456 1'333.086 -0.589
##
## meanCI
## 5'045.879
## 5'306.055
##
## .95
## 6'198.431
##
## kurt
## -0.128
##
## lowest : 2'586.207, 2'937.063, 2'972.973, 3'142.857, 3'178.295
## highest: 6'481.481, 6'614.493, 6'818.530, 6'982.703, 7'333.831