<jupyter_start><jupyter_text>---title: "Exploring Data and Spatial Relationships"output: html_document---<jupyter_code>pre code, pre, code {
  white-space: pre !important;
  overflow-x: auto;
  word-break: keep-all !important;
  word-wrap: initial !important;
}
options(width=800)
options(scipen = 999)<jupyter_output><empty_output><jupyter_text>Learning ObjectivesBy the end of this practical lab you will be able to:* Calculate some descriptive statistics and plot distributions* Perform some basic data cleaning and identify missing values* Explore correlations between variables* Plot maps showing the distribution and spatial patterns of observations  Exploratory data analysisAny type of analysis needs initially some kind of data exploration, which is often very helpful in understanding the underlying patterns of our observations, identify irregularities and selecting the right variables. This step is usually carried out by calculating descriptive statistics, plotting data distributions and other visualization tecqniques, a process also commonly reffered to as _"eyeballing the data"_. Data irregularities can easily thwart our efforts to effectively model urban phenomena, and caution is needed before continuing with any modelling activities.Suppose we decide to do some exploration of broadband connectivity patterns in the London region, such as summary statistics and the distribution of broadband speeds for each neighbourhood (in this case called _Lower Super Output Areas_ [(LSOA)](https://www.ons.gov.uk/methodology/geography/ukgeographies/censusgeographyoutput-area-oa). The raw dataset can be freely downloaded from OFCOM [here](https://www.ofcom.org.uk/research-and-data/multi-sector-research/infrastructure-research/connected-nations-2015/downloads).Firstly we will load some necessary libraries that will help visualize and analyse the data.<jupyter_code>library(classInt)
library(RColorBrewer)
library(ggplot2)
library(ggmap)<jupyter_output><empty_output><jupyter_text>We can then load the data into R using the `read.csv` function.<jupyter_code># Load data
broadband.london <- read.csv("data/bband_london.csv")<jupyter_output><empty_output><jupyter_text>A few common functions used to get some basic information about the data is the `head()` function, which outputs the first 5 observations of the dataset, and `summary()` function which provides some useful descriptive statistics for each variable.<jupyter_code># Preview 
head(broadband.london)

# Get some summary statistics
summary(broadband.london)<jupyter_output><empty_output><jupyter_text>Two common distribution visualisations are the histogram and boxplot. We can use `hist()` to plot the distribution of a variable. As in most functions within R, we can add many parameters to adjust the properties of the histogram, such as the number of density intervals (bins), which we can specify using the `breaks` argument.<jupyter_code># A simple histogram
hist(broadband.london$Median_Downl_Speed, 
     breaks = 15,  # amount of intervals
     density = 20, # optional shading lines
     main = "Distribution of Broadband Download Speeds per LSOA, London", # Graph Title
     xlab = "Median Download Speed (Mbit/s)", # X-Axis Title
     ylab = "Number of LSOAs") # Y-Axis Title<jupyter_output><empty_output><jupyter_text>The next type of graphs that can be useful are _boxplots_. its a type of graph that shows the probability distribution of values along an axis. Boxplots are based on _quartiles_, which divide the dataset into four ranked groups, each representing 25% of the data. Boxplots thus reveal how "spread" values are and highlight any outliers.<jupyter_code># Produce a simple boxplot
boxplot(broadband.london$Median_Downl_Speed,
        main = "Median Download Speeds by LSOA, London",
        ylab = "Median Download Speed (Mbit/s)")<jupyter_output><empty_output><jupyter_text>By default, the "boxed" area of the graph encloses 50% of the values around the median (50% of values, the centre bold line), while the bottom line (hinge) of the box represents the second quartile (25% to 50% of values) and the top at the third (50% to 75% of values). The "whiskers" (the two horizontal lines outside the box) are placed at 1.5 times the Interquartile Range (IQR), which is a value between calculating by difference between the upper (75%) and lower (25%) quartiles, $IQR = Q_3 −  Q_1$. You will also notice that some points are plotted outside the whiskers. If values are above or below the whisker values, i.e. $Q_2 \pm 1.5*IQR$ then these are considered outliers. In this case we have mostly extreme median download speeds on the upper side of the distribution, but they don't seem very suprising. We can get exact information on these values using the `bopxplot.stats()` function.<jupyter_code># Get boxplots statistics
boxplot.stats(broadband.london$Median_Downl_Speed)<jupyter_output><empty_output><jupyter_text>We can use the information provided by `bopxplot.stats()` to add values (labels) to our boxplot, combined with a simple `text()` function after our plot:<jupyter_code># Produce a boxplot with labels
boxplot(broadband.london$Median_Downl_Speed,
        main = "Median Download Speeds by LSOA, London",
        ylab = "Median Download Speed (Mbit/s)")
text(y = boxplot.stats(broadband.london$Median_Downl_Speed)$stats, 
     labels = round(boxplot.stats(broadband.london$Median_Downl_Speed)$stats, 1), 
     x = 1.25) # x controls how far to the right our labels will be plotted<jupyter_output><empty_output><jupyter_text>In R we can easily produce multiple boxplots based the value of another variable, e.g. by type or area. Boxplots in R accept [formulas](https://stat.ethz.ch/R-manual/R-devel/library/stats/html/formula.html), which allow to plot distributions divided by group; in R, the _"by"_ is broadly translated using the `~` symbol, in the fashion of _"value ~ type"_. In this case we will produce separate boxplots for each London subregion and look for any differentiation.<jupyter_code># Boxplots by subregion group
boxplot(Median_Downl_Speed ~ Region, # Formula
        data = broadband.london,
        main = "Median Download Speeds by LSOA, London Subregions",
        xlab = "London Subregion",
        ylab = "Median Download Speed (Mbit/s)")<jupyter_output><empty_output><jupyter_text>The boxplots show the differences in download speed between the different subregions of London. At first glance, differences tend to be subtle but the general trend seems to be that population density may have a key role in broadband speeds. One of the issues here is that grouping areas into subregions may not reflect the actual organization of geographic space. We can increase the granularity and produce a boxplots for each of the 33 London Boroughs.<jupyter_code># Boxplots by Borough 
par(mar = c(11, 4.1, 4.1, 2.1))
boxplot(formula = Median_Downl_Speed ~ Borough, 
        data = broadband.london,
        main = "Median Download Speeds by LSOA, London Boroughs",
        ylab = "Median Download Speed (Mbit/s)", 
        las = 2) # vertical labels
par(mar = c(5.1, 4.1, 4.1, 2.1))<jupyter_output><empty_output><jupyter_text>This graph shows that the differences in the distribution of the Median Download Speeds are much more pronounced at higher scales. There may be dinstinctive local patterns of Internet connectivity; perhaps factors such as population density and income are some of the variables that affect infrastucture development (you can find these data within the same dataset).However, note that increasing granularity offers more insight but diminishes our ability to visualize and easily convey the results; for some, the above graph may already be too difficult to read. Histograms and boxplots are few of the most common tools in visualizing distributions, but there are many other methods that can be used, depending on the nature of data.  CorrelationCorrelation is another concept useful in identifying relationships and exploring urban phenomena. Simply put, correlation is how a variable changes when another variable changes. _"Correlation measures"_ or _"measures of association"_ are the outcomes of statistical methods used to quantify the relationship between two or more variables. Depending on the nature of the variables and the type of relationship various factors / coefficients have been proposed. One of the most common measures is the Pearson's Correlation Coefficient *r*, which measures the the strength and direction of the relationship between two continuous variables. Pearson's r takes values between -1 and 1, with 1 showing a perfect degree of _positive_ association while a -1 a perfect _negative_ respectively. Values around 0 on the other hand suggest that there is a weak relationship between variables.In R, Pearson's Correlation Coefficient can be easily computed using the `cor()` function. For a very basic example, we will use Airbnb data Washington, DC and do some analysis for a snapshot of the offered accommodation. There are various repositories online that can be used to download Airbnb data, such as [insideairbnb.com](http://insideairbnb.com/get-the-data.html) or [tomslee.net](http://tomslee.net/airbnb-data-collection-get-the-data).<jupyter_code># Load data
airbnb.dc <- read.csv("data/airbnb_dc.csv")

# Preview 
head(airbnb.dc)<jupyter_output><empty_output><jupyter_text>Upon inspection, there are some issues with our data that may be problematic in our analysis. For instance, variable _cleaning_fee_ has been identified as a nominal variable (factor) and should be convereted to numerical; generally, numeric values are much easier to handle and model. Empty cells should also be replaced with 0 in this context (although not always!).  There is a dedicated field in statistics dealing with missing data, and this is called _[imputation](https://en.wikipedia.org/wiki/Imputation_(statistics))_. A number of imputation techniques have been proposed, i.e. filling missing values with average values, based on the nearest neighbor or build a model in order to predict values based on other information. Imputation is important since dealing with missing values can be problematic depending on how various models handle them.In this context, missing values in our dataset (R will code those as `NA`) may result in errors depending how various functions handle missing data. Errors also occur because datasets are rarely provided in a format that is ready for analysis. Restructuring data to get to a format that is aligned with the methods of our analysis is an essential and time-wise, often overlooked, step. Remember, data cleaning can take up to 80% - 90% of an analysis. In R, we can retrieve detailed information on data structure with the `str()` function.<jupyter_code># Data structure
str(airbnb.dc)<jupyter_output><empty_output><jupyter_text>We will try fix the _cleaning_fee_ variable to a numeric one:<jupyter_code># Fix cleaning_fee
airbnb.dc$cleaning_fee <- substring(airbnb.dc$cleaning_fee, 2) # crop the first character ($ sign)
airbnb.dc$cleaning_fee[airbnb.dc$cleaning_fee == ""] <- 0 # assign 0 to empty cells
airbnb.dc$cleaning_fee <- as.numeric(airbnb.dc$cleaning_fee)<jupyter_output><empty_output><jupyter_text>Now, suppose we want to perform some analysis on accommodation prices. Let's get some summary statistics for the variables _price_ and _number of beds_.<jupyter_code># Check the distribution of prices and number of beds
summary(airbnb.dc$price)
summary(airbnb.dc$beds)<jupyter_output><empty_output><jupyter_text>Besides descriptive statistics, the `summary` function also gives information about the amount of missing data in our data; in this case 21 NAs are present in the variable _beds_.Remember that outliers can hinder our ability to analyse and visualize distributions. With regards to variable _price_, while the mean is almost $150, the maximum values is $2822. We should avoid plotting extreme outliers as they will obscure details of the distribution; thankfully, we can easily "crop" the visible area of the graph between specified values using the `xlim` and `ylim` for the _x_ and _y_ axis respectively.<jupyter_code># Plot the distribution of price by limiting the x-axis
hist(airbnb.dc$price,
     breaks = max(airbnb.dc$price)/10,  # to get a ~ $10 interval 
     xlim = c(0, 500), # this will limit the x axis values between $0 and $500
     main = "Distribution of Airbnb Prices in Washington, DC, 2015", # Graph Title
     xlab = "Price ($)", # X-Axis Title
     ylab = "Number of Listings") # Y-Axis Title)<jupyter_output><empty_output><jupyter_text>Now, let's calculate the relationship between _price_ and _number of beds_.<jupyter_code># Assign variables for simplicity
x <- airbnb.dc$price
y <- airbnb.dc$beds

# Calculate (Pearson) correlation; note that method = "pearson" is the default, and thus can be omitted
cor(x, y, method = "pearson")<jupyter_output><empty_output><jupyter_text>The above correlation fails because there are NA values in the _beds_ variable. However, we can add the `use = "complete.obs"` argument to exclude observations casewise. Another approach would be to delete observations with NA values from the dataset, altogether but that may remove important information from other variables as well.<jupyter_code>cor(x, y, method = "pearson", use = "complete.obs")<jupyter_output><empty_output><jupyter_text>The r value of ~0.46 suggests that there is a strong positive relationship between price and number of beds, which suggests that accommodation price is associated significantly on the basis of how many people it accommodates. Other factors that may associate with price is, inter alia, housing quality and location. If we were to perform the same analysis in the case of hotel rooms, for example, the number of beds would probably have a very low correlation with price. Since the cost of accommodation is a combination of price (per night) and the cleaning fee, we can also test whether the combined cost is better or worse associated with the number of beds.<jupyter_code># Assign variables
x <- airbnb.dc$price + airbnb.dc$cleaning_fee # add the two variables
y <- airbnb.dc$beds

cor(x, y, method = "pearson", use = "complete.obs")<jupyter_output><empty_output><jupyter_text>Which provides a marginally better result. Still, according to [airbnb](http://rentingyourplace.com/airbnb-101/pricing/to-fee-or-not-to-fee/) the cleaning fee is associated with the duration of stay, while price is per night. In this framework it might not be too sensible to add both of them withough having information on duration of stay. A good practice is to always look at additional information on the variables, usually expressed as _metadata_ accompanying the dataset.With this in mind, lets explore the association between combined price and review score rating.<jupyter_code># Assign variables
x <- airbnb.dc$price + airbnb.dc$cleaning_fee
y <- airbnb.dc$review_scores_rating

cor(x, y, method = "pearson", use = "complete.obs")<jupyter_output><empty_output><jupyter_text>As expected, there is not a strong relationship between review scores and price. The _r_ value is very close to zero, indicating a very weak relationship between the two variables. A "good" accommodation does not mean it is a cheap or expensive one, it depends on multiple factors; to explore such relationships a simple correlation ratio will not suffice, and other techniques, such as regression, should be applied. We will talk more about regression in another section. Note that the strength of a relationship and the significance of a relationship is two different concepts. Usually, correlation is presented with some measure of significance, i.e. the confidence we have that the relationship we identified is not random. We can use `cor.test()` in order to get more information on the test results.<jupyter_code># Assign variables
x <- airbnb.dc$price
y <- airbnb.dc$beds

# Calculate correlation
cor.test(x, y, method = "pearson")<jupyter_output><empty_output><jupyter_text>The test suggest that the correlation is indeed significant at > 95% confidence, as shown by the p-value. Sometimes correlation may be high but not very significant, especially if the number of observations is quite small. It is also important to remember that not all measures of correlation are appropriate for all variables. Pearson correlation is a parametric technique, meaning it will only work correctly if the relationship between variables is linear, both variables are be normally distributed and that values around the regression line are evenly distributed (homoscedasticity). There are many statistical textbooks and [online resources](https://en.wikipedia.org/wiki/Correlation_and_dependence) available that explore these concepts in more detail.The base functions in R includes two more common measures of correlation, Spearman rank correlation and Kendall rank correlation, which are non-parametric, but are designed to calculate strength of association when variables are ordinal. The [UK Index of Multiple Deprivation](https://www.gov.uk/government/statistics/english-indices-of-deprivation-2015), for instance, is rank-based meaning that each area is ranked from worst to best based on various measures of deprivation, so a Pearson correlation between i.e. rank and income would be meaningless.   Spatial AssociationsAnother thing to take into account is that when you compute the correlation coefficient you get a single number; you might however have reason to think that the relationship between variables is not spatially constant. A key assumption correlation statistics is that the observations are independent from one another. This means that the magnitude of a variable at one location should not be correlated with the magnitude at nearby locations. However, this is hardly the case; remember [Tobler's first law of geography](https://en.wikipedia.org/wiki/Tobler%27s_first_law_of_geography). Most urban phenomena appear to have inherent spatial patterns where values are clustered together. This phenomenon is called [_spatial autocorrelation_](https://en.wikipedia.org/wiki/Indicators_of_spatial_association).In our example, the location of accommodation in Washington, DC probably plays an important role in observed prices. For example, central locations might have significantly higher prices than the periphery. The same holds true when we explored broadband connectivity in London earlier.At this stage, we can explore variation by plotting a simple map of price values within Washington, and look or any meaningful spatial patterns. The dataset has two variables, _longitude_ and _latitude_ that we can use to map the accommodation  as points. Variation in prices can be shown using _bins_ of values due to the extreme outliers identified earlier.<jupyter_code># Make bins of prices
breaks <- classIntervals(airbnb.dc$price, 5, style = "kmeans", dataPrecision = 0)

# Make new variable with categories
airbnb.dc$price_cat <- cut(airbnb.dc$price, breaks = breaks$brks, include.lowest = T)

# Get 5 colors based on the YlGnBu ramp
map_colors <- brewer.pal(5, "YlGnBu")<jupyter_output><empty_output><jupyter_text>We can now plot a map of all the accommodation locations colored by price for the Washington, DC area.<jupyter_code># Get the Washington, DC general location
dc.location <- c(-77.11, 38.83, -76.93, 38.99)

# Get backdrop map tiles, various sources can be used here such as google, osm, stamen, etc.
dc.map <- get_map(location=dc.location, source="stamen", maptype="toner", zoom=13, crop=FALSE)

# Make a map with the airbnb points
airbnb.dc.map <- ggmap(dc.map) + 
                 geom_point(aes(x=longitude, y=latitude, color=price_cat), 
                            data=airbnb.dc, alpha = 0.9, size = 0.6) +
                 theme(axis.title=element_blank(), axis.text=element_blank(), 
                 axis.ticks=element_blank(), 
                 panel.border = element_rect(colour="black", fill=NA, size=1.2)) +
                 scale_color_manual(values = map_colors) 

airbnb.dc.map # Plot the map<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.06505520740077589, 'fraction_numerical': 0.007162041181736795, 'mean_word_length': 5.089339794064204, 'pattern_counts': {'":': 0, '<': 108, '<?xml version=': 0, '>': 90, 'https://': 8, 'lorem ipsum': 0, 'www.': 3, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 16, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '22', 'dir': 'jupyter-structured-clean-dedup', 'id': '630424', 'license': 'Apache-2.0', 'n_tokens_mistral': 5471, 'n_tokens_neox': 5118, 'n_words': 2522, 'path': 'notebooks/P12_Spatial_Relationships.ipynb', 'repo_name': 'solartime/UrbanAnalyticsCourse', 'size': '20106'}