PlotUSAStates: Visualization of USA states and their counties

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

View source: R/other_plots.R

Description

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

Usage

1
2
3
4
5
6
7
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)

Arguments

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.

Details

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

Value

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.

Author(s)

Sakis Micheas and Jiaxun Chen

See Also

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

Examples

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
#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)

sppmix documentation built on Jan. 13, 2021, 10:04 p.m.