PlotUSAStates: Visualization of USA states and their counties
In sppmix: Modeling Spatial Poisson and Related Point Processes

Description Usage Arguments Details Value Author(s) See Also Examples

The function plots the requested USA state or county boundaries and additional information if requested or if certain parameters are supplied. We use this function for visualization of geostatistical data, in particular, (Marked) IPPPs.

For examples see

http://faculty.missouri.edu/~micheasa/sppmix/sppmix_all_examples.html#PlotUSAStates

PlotUSAStates(showcounties = FALSE, states = "Missouri",
  showcentroids = TRUE, typecentroid = 0, shownames = FALSE,
  showmarks = FALSE, grayscale = FALSE, open_new_window = FALSE,
  main = "States (true levels)", guidemain = "Level",
  discretelevels = TRUE, levels = 1:3, showplot = TRUE,
  plotlevels = TRUE, marks, pp, surf, boundarycolor = "black",
  namescolor = "black", ppsize = 1)

`showcounties`	Logical to denote that we want a plot of counties. Default is FALSE. Setting this to TRUE will show all the counties for the states passed in the `states` parameter.
`states`	A vector of state names. Set to `NULL` to request all states or `ContinentalUSA_state_names` to show only the continental USA states.
`showcentroids`	Logical requesting to show centroids for each state or county. These centroids are returned in a `ppp` object. The centroid is chosen so that it is always within the state or county boundaries.
`typecentroid`	If `showcentroids=TRUE` we can display either the average of the boundary (`typecentroid=0`) or the "marker point" of the state or county (`typecentroid=1`). For convex states or counties, the latter point is the most south-western point of the state or county.
`shownames`	Logical to display the names of the states for `showcounties=FALSE` or counties for `showcounties=TRUE`.
`showmarks`	Logical to display the mark values given to each state for `showcounties=FALSE` or county for `showcounties=TRUE`.
`grayscale`	Logical to request plots in grayscale.
`open_new_window`	Logical to request plotting in a new graphics window.
`main`	A character string to serve as the main title for the plot.
`guidemain`	A character string to be used as the title for the guide used (legend or colorbar).
`discretelevels`	Logical indicating that the marks are discrete valued.
`levels`	When `discretelevels=TRUE`, the parameter `levels` contains all the possible discrete levels (marks). This is a vector of integers or strings. Default is `1:3`.
`showplot`	Logical requesting to show the plot. Set to FALSE if you want to simply retrieve the centroids of the states or counties, in which case the plot will not be created.
`plotlevels`	Logical requesting that the levels (marks) of each state or county are displayed. If `marks` is not supplied, then for `discretelevels=TRUE` the mark of each state or county is uniformly generated over the values of `levels`, otherwise the marks are uniform in (0,1) (probabilities). If `marks` are given, then they are used to appropriately paint a state or county.
`marks`	A vector of length equal to the number of states or counties requested, containing the mark values for each state or county. A mark is an integer pointing to an element from the vector `levels` for `discretelevels=TRUE`, otherwise a real number.
`pp`	Optionally, a point pattern as an object of type `ppp` to be displayed over the created plot. The window of this point pattern will be used as the window of observation (overrides the window in the `surf` parameter).
`surf`	Optionally, an intensity surface as an object of type `intensity_surface` or an image object of class `im` to be plotted first and then the map will be displayed over this field. Supplying this parameter sets the flag `plotlevels=FALSE` automatically. The window of this intensity surface will be used as the window of observation.
`boundarycolor`	A specific color to use for drawing boundaries. Default is "black". Set to `NULL` if you do not want boundaries drawn.
`namescolor`	A specific color to use for drawing the state or county names when `plotnames=TRUE`. Default is "black".
`ppsize`	Size used in plotting the points. Default is 1.

Note that we use the state and county longitude and latitude boundaries in the USAStatesCounties2016 object.

A list containing the following components:

`PPPcent`	The centroids of the states or counties requested, returned as a marked point pattern.
`PPPMarker`	The marker points of the states or counties requested, returned as a marked point pattern.
`itemnames`	Vector of strings containing all items processed (i.e., either all state names or all county names).
`p`	The created plot, otherwise NULL.

Sakis Micheas and Jiaxun Chen

est_MIPPP_cond_loc, est_mix_damcmc, est_mix_bdmcmc, plot_CompDist, drop_realization, GetBDTable, GetBDCompfit, plotmix_2d, GetBMA, plot_MPP_probs, GetMAPEst

#plot the continental USA with uniformly sampled discrete marks from 10 different levels
ret=PlotUSAStates(states=ContinentalUSA_state_names, levels=1:10, grayscale = FALSE,
 shownames=TRUE, plotlevels =TRUE, discretelevels=TRUE, main="Continental USA (generated levels)")
#now use continuous marks
ret=PlotUSAStates(states=ContinentalUSA_state_names, shownames=FALSE, discretelevels=FALSE,
 main="Continental USA (generated probabilities)", guidemain="Probability", showcentroids = FALSE)
#Fit an IPPP to the California Earthquake data
fitDA=est_mix_damcmc(CAQuakes2014.RichterOver3.0, 8, L = 20000)
#get the surface of Maximum a Posteriori estimates
MAPsurf=GetMAPEst(fitDA)
#plot the states and the earthquake points along with the fitted MAP IPPP intensity surface
ret=PlotUSAStates(states=c('California','Nevada','Arizona'), showcentroids=FALSE,
 shownames=TRUE, main="Earthquakes in CA, 2014", pp=CAQuakes2014.RichterOver3.0, surf=MAPsurf,
 boundarycolor="white", namescolor="white")
#Visualize the Tornado data about MO
#plot the states and the tornado points
ret=PlotUSAStates(states=c('Iowa','Arkansas','Missouri','Illinois','Indiana','Kentucky',
 'Tennessee','Kansas','Nebraska','Texas','Oklahoma','Mississippi','Alabama','Louisiana'),
 showcentroids=FALSE, shownames=TRUE, plotlevels = FALSE, main="Tornadoes about MO, 2011",
 pp=Tornadoes2011MO)
#Visualize aggregate income levels in MO by county using data from the American Community
#Survey (ACS)
#plot in the original scale first; here we pass the marks vector which contains the aggregate
#income values of Missourian counties
ret=PlotUSAStates(showcounties=TRUE, states=c('Missouri'), showcentroids=TRUE, typecentroid=1,
 discretelevels=FALSE, shownames=TRUE, plotlevels=TRUE, marks=MOAggIncomeLevelsPerCounty,
 main="Aggregate Income in MO, 2014", guidemain = "Income level", namescolor="gray",
 boundarycolor="gray")
#plot in the log scale
ret=PlotUSAStates(showcounties=TRUE, states=c('Missouri'), showcentroids=TRUE, typecentroid=1,
 discretelevels=FALSE, shownames=TRUE, plotlevels=TRUE, marks=log(MOAggIncomeLevelsPerCounty),
 main="Aggregate Income in MO, 2014", guidemain = "Income level\n(log scale)", namescolor="gray",
 boundarycolor="gray")
#plot the marker points, county boundaries and names
ret=PlotUSAStates(showcounties=TRUE, states=c('Missouri'), showcentroids=TRUE, typecentroid = 1,
 discretelevels=FALSE, shownames=TRUE, plotlevels=FALSE, marks=log(MOAggIncomeLevelsPerCounty),
 main="Marker points for Missouri counties")
#now plot only the marker points, we treat this as a marked IPPP
ret=PlotUSAStates(showcounties=TRUE, states=c('Missouri'), showcentroids=TRUE, typecentroid = 1,
 discretelevels=FALSE, shownames=FALSE, plotlevels=FALSE, marks=log(MOAggIncomeLevelsPerCounty),
 main="Marker points for Missouri counties", boundarycolor = NULL)
#let us discretize log(income) to 3 levels; low if <=20, average if >20 and <=23, and high if >23
newmarks=rep(0,length(MOAggIncomeLevelsPerCounty))
newmarks[log(MOAggIncomeLevelsPerCounty)<=20]=1
newmarks[log(MOAggIncomeLevelsPerCounty)>20 & log(MOAggIncomeLevelsPerCounty)<=23]=2
newmarks[log(MOAggIncomeLevelsPerCounty)>23]=3
table(newmarks)
levels=c("low","average","high")
ret=PlotUSAStates(showcounties=TRUE, states=c('Missouri'), showcentroids=TRUE, typecentroid=1,
 discretelevels=TRUE, shownames=TRUE, plotlevels=TRUE, main="Aggregate Income in MO, 2014",
 marks=newmarks, levels=levels, guidemain = "Income level", namescolor="gray",
 boundarycolor="gray")
#now fit a marked IPPP model, use the PP of marker points
MPP=ret$PPPMarker
mpp_est <- est_MIPPP_cond_loc(MPP,r=1, hyper=0.2)
plot_MPP_probs(mpp_est)
#now obtain a BDMCMC fit for the ground process this way we can cluster the data
BDMCMCfit <- est_mix_bdmcmc(MPP,m=10,L = 50000)
plot_CompDist(BDMCMCfit)
#use the original output of BDMCMC and apply 10% burnin (default)
BDMCMCfit=drop_realization(BDMCMCfit)
#get the realizations corresponding to the MAP number of components
BDtab=GetBDTable(BDMCMCfit,FALSE)#retrieve frequency table and MAP estimate for
#the number of components
MAPm=BDtab$MAPcomp
BDMCMCfitMAPcomp=GetBDCompfit(BDMCMCfit,MAPm)
BDMCMCfitMAPcompgens=BDMCMCfitMAPcomp$BDgens
MAPsurf=GetMAPEst(BDMCMCfitMAPcompgens)
plotmix_2d(MAPsurf,MPP)+add_title(
 "IPPP intensity surface of MAP estimates (MAP number of components)",
 lambda =MAPsurf$lambda, m=MAPsurf$m, n=MPP$n, L=MAPsurf$L)
plot_ind(BDMCMCfitMAPcompgens)
ret=PlotUSAStates(showcounties=TRUE, states=c('Missouri'),
 showcentroids=TRUE, typecentroid=1, discretelevels=TRUE, shownames=TRUE,
 main="Ground surface of MAP estimates", marks=newmarks, levels=levels,
 guidemain = "Income level", namescolor="gray", boundarycolor="gray",
 pp=MPP, surf=MAPsurf)
#obtain and plot the Bayesian model average; first drop the bad realizations
BDMCMCfit=drop_realization(BDMCMCfit,(BDMCMCfit$Badgen==1))
BMAest=GetBMA(BDMCMCfit)
ret=PlotUSAStates(showcounties=TRUE, states=c('Missouri'),
 showcentroids=TRUE, typecentroid=1, discretelevels=TRUE, shownames=TRUE,
 main="Bayesian model average ground intensity surface", marks=newmarks,
 levels=levels, guidemain = "Income level", namescolor="gray",
 boundarycolor="gray", pp=MPP, surf=BMAest)