date()
## [1] "Wed Apr 16 15:45:05 2014"
This project will investigate spatial patterns in Florida tornado touchdown observations. The extent to which the observed data is explained by clustering due to uneven population density will be explored. A second focus will be how spatial and clustering patterns differ between tornadoes spawned by tropical cyclone (TC) landfalls.
The first step is to load datasets of all tornado touchdowns in the US from 1950-2011, as well as the locations of all TC landfalls over the same timeframe, as well as to load the necessary R packages to perform this study.
loc = "http://moe.met.fsu.edu/~rtruchel/private/all_lf_ondj.txt"
H = read.table(loc, header = TRUE)
download.file("http://myweb.fsu.edu/jelsner/data/tornado2011.zip", "tornado2011.zip",
mode = "wb")
unzip("tornado2011.zip")
library(spatstat)
##
## spatstat 1.36-0 (nickname: 'Intense Scrutiny')
## For an introduction to spatstat, type 'beginner'
library(sp)
library(rgdal)
## rgdal: version: 0.8-16, (SVN revision 498)
## Geospatial Data Abstraction Library extensions to R successfully loaded
## Loaded GDAL runtime: GDAL 1.9.2, released 2012/10/08
## Path to GDAL shared files: /Library/Frameworks/R.framework/Versions/3.0/Resources/library/rgdal/gdal
## Loaded PROJ.4 runtime: Rel. 4.8.0, 6 March 2012, [PJ_VERSION: 480]
## Path to PROJ.4 shared files: /Library/Frameworks/R.framework/Versions/3.0/Resources/library/rgdal/proj
library(rgeos)
## rgeos version: 0.3-3, (SVN revision 437)
## GEOS runtime version: 3.3.3-CAPI-1.7.4
## Polygon checking: TRUE
library(maptools)
## Checking rgeos availability: TRUE
library(ggplot2)
require(maps)
## Loading required package: maps
require(ggmap)
## Loading required package: ggmap
The next step is to prepare the TC data for use as a factor variable in the tornado dataset. In order to do this, all TCs that made landfall in Florida, Alabama, or Mississippi (due to tornadoes predominantly being found in outer TC bands) between 1950-2011 are selected. Dates that are within plus or minus two days from a landfall are considered TC days, and any Florida tornadoes observed in these ranges are considered to be TC-related tornadoes.
min = which.max(H$X1900 > 1949)
max = which.min(H$X1900 < 2012)
H1 = H[(min + 1):(max - 1), ]
mmm = c(0, 31, 59, 90, 120, 151, 181, 212, 243, 273, 304, 334)
H1$Date = (mmm[H1$X09] + H1$X09.1) + 365 * ((H1$X1900) - 1950)
H2 = H1[1, ]
cd = dim(H1)
cd1 = cd[1]
cd2 = cd[2]
i1 = 1
for (i in 1:cd1) {
if (H1$X01.1[i] > 2) {
if (H1$X01.1[i] < 9) {
H2[i1, ] = H1[i, ]
i1 = i1 + 1
}
}
}
cd3 = dim(H2)
cd4 = cd3[1]
lfprox = matrix(0, 22630, 1)
for (j in 1:cd4) {
j1 = H2$Date[j]
j2 = ((j1 - 2):(j1 + 2))
lfprox[j2, 1] = 1
}
This TC or non-TC classification is now added to the table of tornado data as a new column, by indexing the date information in the tornado dataset against the TC and non-TC days prepared in the previous step.
Torn = readOGR(dsn = ".", layer = "tornado", stringsAsFactors = FALSE)
## OGR data source with driver: ESRI Shapefile
## Source: ".", layer: "tornado"
## with 56221 features and 28 fields
## Feature type: wkbPoint with 2 dimensions
Torn1 = subset(Torn, STATE == "FL")
Torn1$Date = (mmm[Torn1$MONTH] + Torn1$DAY) + 365 * ((Torn1$YEAR) - 1950)
Torn1$TC = 0
cd5 = dim(Torn1)
cd6 = cd5[1]
for (k in 1:cd6) {
Torn1$TC[k] = lfprox[Torn1$Date[k]]
}
The table of Florida tornadoes can now be divided into four different subsets for further study: All Non-TC related tornadoes, all TC tornadoes, weak (EF0/1) non-TC tornadoes, and strong (EF2+) non-TC tornadoes.
TornNTC = subset(Torn1, TC == 0)
TornTC = subset(Torn1, TC == 1)
TornWeak = subset(Torn1, FSCALE <= 1)
TornStrong = subset(Torn1, FSCALE > 1)
TornNTC$EF = factor(TornNTC$FSCALE)
TornTC$EF = factor(TornTC$FSCALE)
TornWeak$EF = factor(TornWeak$FSCALE)
TornStrong$EF = factor(TornStrong$FSCALE)
TornTC0 = subset(Torn1, FSCALE == 0 & TC == 1)
TornTC1 = subset(Torn1, FSCALE == 1 & TC == 1)
TornTC2 = subset(Torn1, FSCALE == 2 & TC == 1)
TornTC3 = subset(Torn1, FSCALE == 3 & TC == 1)
TornTC4 = subset(Torn1, FSCALE == 4 & TC == 1)
## Warning: no non-missing arguments to min; returning Inf
## Warning: no non-missing arguments to max; returning -Inf
## Warning: no non-missing arguments to min; returning Inf
## Warning: no non-missing arguments to max; returning -Inf
TornTC5 = subset(Torn1, FSCALE == 5 & TC == 1)
## Warning: no non-missing arguments to min; returning Inf
## Warning: no non-missing arguments to max; returning -Inf
## Warning: no non-missing arguments to min; returning Inf
## Warning: no non-missing arguments to max; returning -Inf
TornNTC0 = subset(Torn1, FSCALE == 0 & TC == 0)
TornNTC1 = subset(Torn1, FSCALE == 1 & TC == 0)
TornNTC2 = subset(Torn1, FSCALE == 2 & TC == 0)
TornNTC3 = subset(Torn1, FSCALE == 3 & TC == 0)
TornNTC4 = subset(Torn1, FSCALE == 4 & TC == 0)
TornNTC5 = subset(Torn1, FSCALE == 5 & TC == 0)
## Warning: no non-missing arguments to min; returning Inf
## Warning: no non-missing arguments to max; returning -Inf
## Warning: no non-missing arguments to min; returning Inf
## Warning: no non-missing arguments to max; returning -Inf
EFCount = matrix(0, 2, 6)
EFCount[1, 1] = (dim(TornTC0)[1])
EFCount[1, 2] = (dim(TornTC1)[1])
EFCount[1, 3] = (dim(TornTC2)[1])
EFCount[1, 4] = (dim(TornTC3)[1])
EFCount[1, 5] = (dim(TornTC4)[1])
EFCount[1, 6] = (dim(TornTC5)[1])
EFCount[2, 1] = (dim(TornNTC0)[1])
EFCount[2, 2] = (dim(TornNTC1)[1])
EFCount[2, 3] = (dim(TornNTC2)[1])
EFCount[2, 4] = (dim(TornNTC3)[1])
EFCount[2, 5] = (dim(TornNTC4)[1])
EFCount[2, 6] = (dim(TornNTC5)[1])
EFCount
## [,1] [,2] [,3] [,4] [,5] [,6]
## [1,] 215 79 39 2 0 0
## [2,] 1542 758 276 33 2 0
barplot(EFCount, main = "Florida tornado counts, 1950-2011, \nby TC (blue) and non-TC (red) origins",
ylab = "Tornado Count", col = c("darkblue", "red"), names.arg = c("EF0",
"EF1", "EF2", "EF3", "EF4", "EF5"), beside = TRUE)
The resulting bar plot shows that only around 11% of Florida tornadoes are TC-related, with a roughly equal distribution between TC and non-TC tornadoes in terms of EF rating. Florida tornadoes tend to be weak, with only 12% rated as EF2 or higher. As weak and short-lived tornadoes are more likely to not be observed in unpopulated areas, this may contribute to more population clustering in Florida tornadoes.
The next step will be to spatially analyze the data, first by looking at regional count and density data.
counties = maps::map("county", "florida", fill = TRUE, plot = FALSE)
ll = "+proj=longlat +datum=WGS84"
IDs = sapply(strsplit(counties$names, ","), function(x) x[2])
counties.spll = maptools::map2SpatialPolygons(counties, IDs = IDs, proj4string = CRS(ll))
counties.spT = sp::spTransform(counties.spll, CRS(proj4string(TornNTC)))
centers.spT = rgeos::gCentroid(counties.spT, byid = TRUE)
centers.spll = sp::spTransform(centers.spT, CRS(ll))
centers.dfll = as.data.frame(coordinates(centers.spll))
names(centers.dfll) = c("long", "lat")
all(rownames(centers.dfll) == names(counties.spT))
## [1] TRUE
ct = over(counties.spT, TornNTC, returnList = TRUE)
nT = sapply(ct, nrow)
county = rownames(centers.dfll)
area = rgeos::gArea(counties.spT, byid = TRUE)
nTor.df = data.frame(state = "FL", county = county, nT = nT, area = area, rate = nT/area *
10^6 * 100, stringsAsFactors = FALSE)
spdf = SpatialPolygonsDataFrame(counties.spT, nTor.df)
ct1 = over(counties.spT, TornTC, returnList = TRUE)
nT1 = sapply(ct1, nrow)
nTor1.df = data.frame(state = "FL", county = county, nT1 = nT1, area = area,
rate = nT1/area * 10^6 * 100, stringsAsFactors = FALSE)
spdf1 = SpatialPolygonsDataFrame(counties.spT, nTor1.df)
ct2 = over(counties.spT, TornWeak, returnList = TRUE)
nT2 = sapply(ct2, nrow)
nTor2.df = data.frame(state = "FL", county = county, nT2 = nT2, area = area,
rate = nT2/area * 10^6 * 100, stringsAsFactors = FALSE)
spdf2 = SpatialPolygonsDataFrame(counties.spT, nTor2.df)
ct3 = over(counties.spT, TornStrong, returnList = TRUE)
nT3 = sapply(ct3, nrow)
nTor3.df = data.frame(state = "FL", county = county, nT3 = nT3, area = area,
rate = nT3/area * 10^6 * 100, stringsAsFactors = FALSE)
spdf3 = SpatialPolygonsDataFrame(counties.spT, nTor3.df)
spplot(spdf, "nT", main = "Florida non-TC-related tornado touchdowns by county")
spplot(spdf1, "nT1", main = "Florida TC-Related Tornado Touchdowns by County")
spplot(spdf2, "nT2", main = "Florida Weak Tornado Touchdowns by County")
spplot(spdf3, "nT3", main = "Florida Strong Tornado Touchdowns by County")
fl.df = map_data("state", "florida")
fl.m = cbind(fl.df$long, fl.df$lat)
fl.P = Polygon(fl.m)
fl.PS = Polygons(list(fl.P), 1)
fl.SPS = SpatialPolygons(list(fl.PS))
proj4string(fl.SPS) = CRS("+proj=longlat +ellps=WGS84")
fl.SPS = spTransform(fl.SPS, CRS(proj4string(TornNTC)))
gIsValid(fl.SPS)
## [1] TRUE
xy = fortify(fl.SPS)
W = owin(poly = list(x = rev(xy$long[-1]), y = rev(xy$lat[-1])))
TornNTC.ppp = as.ppp(TornNTC["EF"])
TornNTCFL.ppp = TornNTC.ppp[W]
unitname(TornNTCFL.ppp) = c("meter", "meters")
fl.SPS = spTransform(fl.SPS, CRS(proj4string(TornTC)))
gIsValid(fl.SPS)
## [1] TRUE
xy = fortify(fl.SPS)
TornTC.ppp = as.ppp(TornTC["EF"])
TornTCFL.ppp = TornTC.ppp[W]
unitname(TornTCFL.ppp) = c("meter", "meters")
summary(TornTCFL.ppp)
## Marked planar point pattern: 265 points
## Average intensity 1.79e-09 points per square meter
##
## *Pattern contains duplicated points*
##
## Coordinates are given to 2 decimal places
## i.e. rounded to the nearest multiple of 0.01 meters
##
## Multitype:
## frequency proportion intensity
## -9 9 0.034 6.08e-11
## 0 169 0.638 1.14e-09
## 1 57 0.215 3.85e-10
## 2 30 0.113 2.03e-10
## 3 0 0.000 0.00e+00
##
## Window: polygonal boundary
## single connected closed polygon with 871 vertices
## enclosing rectangle: [801878, 1615764] x [-1418362, -799376] meters
## Window area = 1.48021e+11 square meters
## Unit of length: 1 meter
den = density(TornNTCFL.ppp[TornNTCFL.ppp$marks == 1], dimyx = c(256, 256),
adjust = 0.8)
plot(den, main = "Spatial Density of non-TC-related \nTornado Touchdowns in Florida, 1950-2011")
points(TornNTCFL.ppp, pch = 16, cex = 0.75)
den1 = density(TornTCFL.ppp[TornTCFL.ppp$marks == 1], dimyx = c(256, 256), adjust = 0.8)
plot(den1, main = "Spatial Density of Recorded TC-Related \nTornado Touchdowns in Florida, 1950-2011")
points(TornTCFL.ppp, pch = 16, cex = 0.75)
County tornado counts are broken out 4 ways: all non-TC related tornadoes, all TC-related tornadoes, weak tornadoes, and strong tornadoes. In general, the non-TC related tornadoes and weak tornadoes show very few differences between them, and while there appear to be local maxes of strong tornadoes along the I-4 corridor, in Broward and Palm Beach Counties, and in the Western Panhandle, the strong tornado counts are more erratic as they are rare events. For this reason, the rest of the project will focus on the differences between non-TC and TC tornadoes for the purpose of clarity.
The spatial density plots of non-TC and TC tornadoes do reveal that there are notable differences between the patterns of occurance. Non-TC tornado densities are highest along the I-4 corridor, especially in the Tampa Bay region. Additional maxima are located in the Pensacola area as well as coastal southeastern Florida. In contrast, TC tornado density is highest along the East Coast, particularly south of Daytona Beach. Density is elevated along I-4 , but there are notable minima in Southwest Florida and the Big Bend. There is also a strong maximum north of Apalachicola. These patterns are consistent with common TC landfall tracks across Southeast Florida and the northern Gulf Coast, which place these regions in the right front moving quadrant of the incoming TC, which is favored for tornado development.
plot(envelope(TornNTCFL.ppp, fun = Gest, nsim = 99), main = "Nearest Neighbor Distance Function for non-TC-related \n1950-2011 Florida Tornadoes, with Simulated Clustering Envelope")
## Generating 99 simulations of CSR ...
## 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, 86, 87, 88, 89, 90,
## 91, 92, 93, 94, 95, 96, 97, 98, 99.
##
## Done.
## lty col key label
## obs 1 1 obs hat(G)[obs](r)
## theo 2 2 theo G[theo](r)
## hi 1 8 hi hat(G)[hi](r)
## lo 1 8 lo hat(G)[lo](r)
## meaning
## obs observed value of G(r) for data pattern
## theo theoretical value of G(r) for CSR
## hi upper pointwise envelope of G(r) from simulations
## lo lower pointwise envelope of G(r) from simulations
plot(envelope(TornTCFL.ppp, fun = Gest, nsim = 99), main = "Nearest Neighbor Distance Function for TC-related \n1950-2011 Florida Tornadoes, with Simulated Clustering Envelope")
## Generating 99 simulations of CSR ...
## 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, 86, 87, 88, 89, 90,
## 91, 92, 93, 94, 95, 96, 97, 98, 99.
##
## Done.
## lty col key label
## obs 1 1 obs hat(G)[obs](r)
## theo 2 2 theo G[theo](r)
## hi 1 8 hi hat(G)[hi](r)
## lo 1 8 lo hat(G)[lo](r)
## meaning
## obs observed value of G(r) for data pattern
## theo theoretical value of G(r) for CSR
## hi upper pointwise envelope of G(r) from simulations
## lo lower pointwise envelope of G(r) from simulations
Performing a nearest neighbor distance function analysis of the Florida tornado data reveals some major differences between the TC and non-TC cases. As the simulation envelope does not overlap with the observed nearest neighbor distribution of tornadoes for the non-TC-related tornadoes, there is ample evidence to support the existence of clustering in these cases. This means that the location of tornado touchdowns shows evidence of observation bias and a partiality towards populated areas. However, for the TC-related tornadoes, the simulation envelope overlaps the observed distribution for approximately half the radius domain. This means that while there is some evidence for clustering, especially on highly localized scales, in general, TC-related tornadoes tend to be more comprehensively observed. This theory will be further tested in the next step of the analysis.
download.file("http://myweb.fsu.edu/jelsner/data/ci08au12.zip", "ci08au12.zip",
mode = "wb")
unzip("ci08au12.zip")
US = readOGR(dsn = ".", layer = "ci08au12", p4s = "+proj=longlat +ellps=WGS84")
## OGR data source with driver: ESRI Shapefile
## Source: ".", layer: "ci08au12"
## with 43603 features and 23 fields
## Feature type: wkbPoint with 2 dimensions
UST = spTransform(US, CRS(proj4string(Torn)))
USt = subset(UST, UST$POP_1990 > 2500)
US.ppp = as.ppp(USt)
W = as.owin(TornNTCFL.ppp)
towns.ppp = unmark(US.ppp[W])
Zc = distmap(towns.ppp)
rhat = rhohat(TornNTCFL.ppp[TornNTCFL.ppp$marks == 1], Zc, adjust = 2, smoother = "kernel",
method = "transform")
## Warning: data contain duplicated points
m2km = 0.001
km2 = 1e+06
dist = rhat$Zc * m2km
rho = rhat$rho * km2
hi = rhat$hi * km2
lo = rhat$lo * km2
Rho.df = data.frame(dist = dist, rho = rho, hi = hi, lo = lo)
W1 = as.owin(TornTCFL.ppp)
towns1.ppp = unmark(US.ppp[W1])
Zc1 = distmap(towns1.ppp)
rhat1 = rhohat(TornTCFL.ppp[TornTCFL.ppp$marks == 1], Zc1, adjust = 2, smoother = "kernel",
method = "transform")
## Warning: data contain duplicated points
dist1 = rhat1$Zc * m2km
rho1 = rhat1$rho * km2
hi1 = rhat1$hi * km2
lo1 = rhat1$lo * km2
Rho1.df = data.frame(dist1 = dist1, rho1 = rho1, hi1 = hi1, lo1 = lo1)
ggplot(Rho.df) + geom_ribbon(aes(x = dist, ymin = lo, ymax = hi), alpha = 0.3) +
geom_line(aes(x = dist, y = rho), col = "black") + ylab(expression(paste(hat(rho),
"(Tornado Reports/km", {
}^2, ")"))) + xlab("Distance from Nearest Town Center (km)") + ggtitle("Spatial Density of Florida non-TC Tornado Reports \nas a Function of Distance from Town Centers") +
theme_bw()
ggplot(Rho1.df) + geom_ribbon(aes(x = dist1, ymin = lo1, ymax = hi1), alpha = 0.3) +
geom_line(aes(x = dist1, y = rho1), col = "black") + ylab(expression(paste(hat(rho1),
"(Tornado Reports/km", {
}^2, ")"))) + xlab("Distance from Nearest Town Center (km)") + ggtitle("Spatial Density of Florida TC-related Tornado Reports \nas a Function of Distance from Town Centers") +
theme_bw()
Comparing the spatial density of TC and non-TC tornadoes as a function of distance from town centers, there continues to be a clear divergence between the spatial observation patterns of the two datasets. Non-TC tornado observation frequency declines to less than 1/7th of peak value at a distance of 15 km from the nearest town center before rebounding. In contrast, TC tornado density declines to a little less than half of peak values by 15km and flattens out thereafter. This also supports that TC tornado observations are less biased by populated areas overall. This will be explicitly checked in the final section.
marks(towns.ppp) = NULL
marks(TornNTCFL.ppp) = NULL
TornTown.ppp = superimpose(towns.ppp, TornNTCFL.ppp)
## Warning: data contain duplicated points
marks(TornTown.ppp) = as.factor(c(rep("city", towns.ppp$n), rep("tornado", TornNTCFL.ppp$n)))
marks(towns1.ppp) = NULL
marks(TornTCFL.ppp) = NULL
TornTown1.ppp = superimpose(towns1.ppp, TornTCFL.ppp)
## Warning: data contain duplicated points
marks(TornTown1.ppp) = as.factor(c(rep("city", towns1.ppp$n), rep("tornado",
TornTCFL.ppp$n)))
plot(envelope(TornTown.ppp, correction = "border", fun = Kcross, nsim = 99),
main = "Kcross Function with Simulation Envelope, \nFlorida Towns (Pop. >2500) and non-TC-related Tornadoes",
xlim = c(0, 90000))
## Generating 99 simulations of CSR ...
## 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, 86, 87, 88, 89, 90,
## 91, 92, 93, 94, 95, 96, 97, 98, 99.
##
## Done.
## lty col key label
## obs 1 1 obs {hat(K)[list(city,tornado)]^{obs}}(r)
## theo 2 2 theo {K[list(city,tornado)]^{theo}}(r)
## hi 1 8 hi {hat(K)[list(city,tornado)]^{hi}}(r)
## lo 1 8 lo {hat(K)[list(city,tornado)]^{lo}}(r)
## meaning
## obs observed value of Kcross["city", "tornado"](r) for data pattern
## theo theoretical value of Kcross["city", "tornado"](r) for CSR
## hi upper pointwise envelope of Kcross["city", "tornado"](r) from simulations
## lo lower pointwise envelope of Kcross["city", "tornado"](r) from simulations
plot(envelope(TornTown1.ppp, correction = "border", fun = Kcross, nsim = 99),
main = "Kcross Function with Simulation Envelope, \nFlorida Towns (Pop. >2500) and TC-related Tornadoes",
xlim = c(0, 90000))
## Generating 99 simulations of CSR ...
## 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, 86, 87, 88, 89, 90,
## 91, 92, 93, 94, 95, 96, 97, 98, 99.
##
## Done.
## lty col key label
## obs 1 1 obs {hat(K)[list(city,tornado)]^{obs}}(r)
## theo 2 2 theo {K[list(city,tornado)]^{theo}}(r)
## hi 1 8 hi {hat(K)[list(city,tornado)]^{hi}}(r)
## lo 1 8 lo {hat(K)[list(city,tornado)]^{lo}}(r)
## meaning
## obs observed value of Kcross["city", "tornado"](r) for data pattern
## theo theoretical value of Kcross["city", "tornado"](r) for CSR
## hi upper pointwise envelope of Kcross["city", "tornado"](r) from simulations
## lo lower pointwise envelope of Kcross["city", "tornado"](r) from simulations
Computing an envelope of simulated K-Cross functions and comparing with the observed K-Cross function confirms that the populations of non-TC and TC tornadoes are distinct. As the simulation envelope and the observed K-Cross function mostly do not overlap (first above, then below) for the non-TC tornadoes, it appears that much of the clustering seen in non-TC tornado touchdowns in Florida is in fact explained by city locations and the resulting observation bias. However, for TC tornadoes, the envelope and observed K-Cross functions mostly overlap. This indicates that while population clustering may have some explanatory power at particular observational radii, for the most part spatial patterns in TC tornado touchdowns are not explained simply by city locations in Florida.
Conclusion: Across a number of statistical tests of clustering, non-TC and TC tornadoes observed since 1950 in Florida do exhibit different behaviors. Overall, observations non-TC tornadoes appear to be strongly clustered in a way that is explained by population patterns, indicating relatively more events may have been missed. TC-related tornadoes do not exhibit such tendencies, so fewer of these type of tornadoes may be missing from the historical record of touchdowns. This may possibly be because Florida's coastal regions are favored for TC-related tornadoes, where without respect to population density it may be easier to observe tornado damage patterns or spot approaching funnel clouds. In any case, the observed differences between these two “flavors” or tornado in Florida and other coastal states is an interesting question for further study.