Project Specifications

Texas Realty Insights aims to analyze real estate market trends in the state of Texas by analyzing historical real estate sales data. The goal is to provide statistical and visual insights to support strategic decisions regarding sales and the optimization of real estate advertisements.

The objectives of the project are as follows:

  • Identify and interpret historical trends in real estate sales in Texas.
  • Evaluate the efficacy of marketing strategies for real estate listings.
  • Provide a graphical visualization of the data that shows the distribution of prices and sales across cities, months, and years.

The dataset contains the following variables:

  • city: reference city
  • year: reference year
  • month: reference month
  • sales: total number of sales
  • volume: total sales value (in millions of dollars)
  • median_price: median sales price (in dollars)
  • listings: total number of active advertisements
  • months_inventory: the amount of time required to sell all current ads, expressed in months.

1 Analysis of variables

First of all, we want to improve the presentation of the tables by defining a dedicated function using kable.

# ----------------------- TABLE FORMATTING FUNCTION -----------------------

show_table <- function(data, digits = 2, scrollable = FALSE){
  # data: input data frame to display
  # digits: number of decimal places to show in numeric columns
  # scrollable: logical, whether to enable horizontal scrolling
  
  numeric_cols <- sapply(data, is.numeric)
  data[numeric_cols] <- lapply(data[numeric_cols], round, digits = digits)
  
  table <- kable(data) %>%
    kable_styling(
      bootstrap_options = c("striped", "hover", "condensed", "bordered"),
      full_width = FALSE,
      position = "center"
    )
  
  if (scrollable == TRUE) {
    table <- table %>%
      scroll_box(width = "100%", height = "auto")
  }
  
  return(table)
}

# -------------------------------------------------------------------------

Now we load the CSV dataset from the given URL, inspect its dimensions, verify the presence of missing values, and display the first ten records.

# Load dataset
dataset_url <- "https://drive.google.com/uc?export=download&id=1O4If8876MTwstkrZX0BqpQ_BxcsIMEko"
dataset <- read.csv(dataset_url)

# Inspect the dimensions of dataset
records <- dim(dataset)[1]
variables <- dim(dataset)[2]

# Verify whether the dataset contains any missing values
cat("The dataset", 
    ifelse(anyNA(dataset), "contains", "does not contain"), 
    "missing values.\n")

# Display the first 10 records
show_table(head(dataset, 10))
## The dataset does not contain missing values.
city year month sales volume median_price listings months_inventory
Beaumont 2010 1 83 14.16 163800 1533 9.5
Beaumont 2010 2 108 17.69 138200 1586 10.0
Beaumont 2010 3 182 28.70 122400 1689 10.6
Beaumont 2010 4 200 26.82 123200 1708 10.6
Beaumont 2010 5 202 28.83 123100 1771 10.9
Beaumont 2010 6 189 27.22 122800 1803 11.1
Beaumont 2010 7 164 22.71 124300 1857 11.7
Beaumont 2010 8 174 25.24 136800 1830 11.6
Beaumont 2010 9 124 17.23 121100 1829 11.7
Beaumont 2010 10 150 23.90 138500 1779 11.5

The dataset contains 8 variables and 240 records. Let’s show the structure of the dataset and the types of variables.

# Show the structure of the dataset
str(dataset)
## 'data.frame':    240 obs. of  8 variables:
##  $ city            : chr  "Beaumont" "Beaumont" "Beaumont" "Beaumont" ...
##  $ year            : int  2010 2010 2010 2010 2010 2010 2010 2010 2010 2010 ...
##  $ month           : int  1 2 3 4 5 6 7 8 9 10 ...
##  $ sales           : int  83 108 182 200 202 189 164 174 124 150 ...
##  $ volume          : num  14.2 17.7 28.7 26.8 28.8 ...
##  $ median_price    : num  163800 138200 122400 123200 123100 ...
##  $ listings        : int  1533 1586 1689 1708 1771 1803 1857 1830 1829 1779 ...
##  $ months_inventory: num  9.5 10 10.6 10.6 10.9 11.1 11.7 11.6 11.7 11.5 ...

We can identify the types of statistical variables in the dataset as follows:

  • city: categorical nominal variable.
  • year: categorical ordinal variable or discrete quantitative variable on an interval scale, depending on the analysis.
  • month: categorical ordinal variable.
  • sales: discrete quantitative variable on a ratio scale.
  • volume: continuous quantitative variable on a ratio scale.
  • median_price: continuous quantitative variable on a ratio scale.
  • listings: discrete quantitative variable on a ratio scale.
  • months_inventory: continuous quantitative variable on a ratio scale.

We can carry out different types of analysis depending on the variables in the dataset:

  • Geographical analysis of sales, volume, median price and listings across cities, in order to rank cities according to Texas real estate market performance.
  • Time-series analysis of sales, volume, median price and listings, to study how these variables evolve over time.
  • Monthly and seasonal analysis to determine whether certain months or seasons perform better than others.

Before starting our analysis, let’s attach the dataset so that the variables are accessible in the workspace.

attach(dataset)

Now we want to handle the time-based variables more effectively. To do so, we will add to the dataset a date column in day-month-year format, setting the day to the first day of each month, so that we can later perform time series analysis.

# Add the date column to the dataset
dataset$date <- as.Date(paste("01", month, year, sep = "-"), format = "%d-%m-%Y")

Now we want to convert the months (labeled as \(1\), \(2\), \(3\), etc.) into text month names for clearer visualization.

dataset$month <- factor(
  dataset$month,
  levels = 1:12,
  labels = c("Jan", "Feb", "Mar", "Apr", "May", "Jun",
             "Jul", "Aug", "Sep", "Oct", "Nov", "Dec")
)

The dataset is now ready for analysis. Let’s inspect its final structure.

# Show the structure of the dataset
str(dataset)
## 'data.frame':    240 obs. of  9 variables:
##  $ city            : chr  "Beaumont" "Beaumont" "Beaumont" "Beaumont" ...
##  $ year            : int  2010 2010 2010 2010 2010 2010 2010 2010 2010 2010 ...
##  $ month           : Factor w/ 12 levels "Jan","Feb","Mar",..: 1 2 3 4 5 6 7 8 9 10 ...
##  $ sales           : int  83 108 182 200 202 189 164 174 124 150 ...
##  $ volume          : num  14.2 17.7 28.7 26.8 28.8 ...
##  $ median_price    : num  163800 138200 122400 123200 123100 ...
##  $ listings        : int  1533 1586 1689 1708 1771 1803 1857 1830 1829 1779 ...
##  $ months_inventory: num  9.5 10 10.6 10.6 10.9 11.1 11.7 11.6 11.7 11.5 ...
##  $ date            : Date, format: "2010-01-01" "2010-02-01" ...

2 Measures of central tendency, variability and shape

In this section we aim to compute the measures of central tendency, dispersion and shape for the quantitative variables (sales, volume, median_price, listings and months_inventory).

For the other variables (city, year, month), which can all be considered categorical variables, we will show the tables of absolute and relative frequencies of their categories.

2.1 Statistical measures for quantitative variables

Let’s create a function that displays the statistical measures all in one. In particular, we want to evaluate:

  • Measures of central tendency: mean, minimum, first quartile, median, third quartile, maximum;
  • Measures of variability: variance, standard deviation, range, interquartile range, coefficient of variation;
  • Measures of shape: Fisher’s skewness and excess kurtosis.
# ---------------------- STATISTICAL SUMMARY FUNCTION ----------------------

stat_summary <- function(df, variables){
  # df: data frame containing the data to analyze
  # variables: names of the columns to summarize
  
  results <- list()
  
  for(v in variables){
    x <- df[[v]]
    
    if(is.numeric(x)){
      mu <- mean(x, na.rm = TRUE)
      sigma <- sd(x, na.rm = TRUE)
      min_val <- min(x, na.rm = TRUE)
      max_val <- max(x, na.rm = TRUE)
      
      results[[v]] <- c(
        # Measures of central tendency
        mean = mu,
        min = min_val,
        q1 = as.numeric(quantile(x, probs = 0.25, na.rm = TRUE)),
        median = median(x, na.rm = TRUE),
        q3 = as.numeric(quantile(x, probs = 0.75, na.rm = TRUE)),
        max = max_val,
        
        # Measures of variability
        variance = var(x, na.rm = TRUE),
        std_dev = sigma,
        range = max_val - min_val,
        interquartile_range = IQR(x, na.rm = TRUE),
        cv = ifelse(mu == 0, NA, sigma/mu),
  
        # Measures of shape
        skewness_index = skewness(x, na.rm = TRUE),
        kurtosis_index = kurtosis(x, na.rm = TRUE) - 3
      )
    }
  }
  
  # Convert the list of results into a data frame
  stats <- as.data.frame(results)
  
  # Assign names to the statistics
  rownames(stats) <- c(
  "Mean value", "Minimum", "First quartile (Q1)", "Median value", 
  "Third quartile (Q3)", "Maximum", "Variance", "Standard deviation", 
  "Range", "Interquartile range", "Coefficient of variation",
  "Fisher's Skewness", "Excess Kurtosis"
  )
  
  return(stats)
}

# --------------------------------------------------------------------------
# Compute statistical measures for the selected variables
stats <- stat_summary(dataset, c("sales", "volume", "median_price", "listings", "months_inventory"))

# Display the statistical summary rounded with 2 decimal
show_table(stats)
sales volume median_price listings months_inventory
Mean value 192.29 31.01 132665.42 1738.02 9.19
Minimum 79.00 8.17 73800.00 743.00 3.40
First quartile (Q1) 127.00 17.66 117300.00 1026.50 7.80
Median value 175.50 27.06 134500.00 1618.50 8.95
Third quartile (Q3) 247.00 40.89 150050.00 2056.00 10.95
Maximum 423.00 83.55 180000.00 3296.00 14.90
Variance 6344.30 277.27 513572983.09 566568.97 5.31
Standard deviation 79.65 16.65 22662.15 752.71 2.30
Range 344.00 75.38 106200.00 2553.00 11.50
Interquartile range 120.00 23.23 32750.00 1029.50 3.15
Coefficient of variation 0.41 0.54 0.17 0.43 0.25
Fisher’s Skewness 0.72 0.88 -0.36 0.65 0.04
Excess Kurtosis -0.31 0.18 -0.62 -0.79 -0.17

As we can see, the quantitative variables contain information on very different scales and orders of magnitude. Let’s summarize some important observations:

  • The coefficients of variation suggest substantial dispersion for some variables, especially sales, volume and listings, while the median_price is more stable.
  • The skewness values, evaluated using Fisher’s index, indicate that all variables are not symmetric, with the exception of months_inventory, which has a skewness near to zero.
  • The kurtosis values suggest deviations from Gaussian shape, with distributions that are often platykurtic.

2.2 Distribution of frequencies for categorical variables

Now let’s proceed creating a function that displays the distribution of frequencies for the qualitative variables.

# -------------------- FREQUENCY DISTRIBUTION FUNCTION --------------------

freq_distribution <- function(df, variable){
  # df: dataframe containing the data to analyze
  # variable: name of the column to evaluate
  
  x <- df[[variable]]
  N <- length(x)
  ni <- table(x)   # absolute frequencies
  fi <- ni/N       # relative frequencies
  
  freq_summary <- data.frame(
    labels = names(ni),
    abs_freq = as.vector(ni),
    rel_freq = as.vector(fi)
  )
  
  colnames(freq_summary) <- c(variable, "Frequency", "Relative frequency")
    
  return(freq_summary)
}

# -------------------------------------------------------------------------

Now we can visualize the frequency tables for the qualitative variables.

# Frequency distribution of city
freqs_city <- freq_distribution(dataset, "city")
show_table(freqs_city)
city Frequency Relative frequency
Beaumont 60 0.25
Bryan-College Station 60 0.25
Tyler 60 0.25
Wichita Falls 60 0.25 
# Frequency distribution of year
freqs_year <- freq_distribution(dataset, "year")
show_table(freqs_year)
year Frequency Relative frequency
2010 48 0.2
2011 48 0.2
2012 48 0.2
2013 48 0.2
2014 48 0.2 
# Frequency distribution of month
freqs_month <- freq_distribution(dataset, "month")
show_table(freqs_month)
month Frequency Relative frequency
Jan 20 0.08
Feb 20 0.08
Mar 20 0.08
Apr 20 0.08
May 20 0.08
Jun 20 0.08
Jul 20 0.08
Aug 20 0.08
Sep 20 0.08
Oct 20 0.08
Nov 20 0.08
Dec 20 0.08

We can observe that the frequency distribution is uniform within each qualitative variable, as all categories occur with the same frequency in the dataset.

3 Variables with the highest variability and asymmetry

Since the quantitative variables have very different scales, an useful estimator for measuring variability is the coefficient of variation: it is calculated as the standard deviation divided by the mean (\(CV = \frac{\sigma}{\bar{x}}\)), and it measures the relative dispersion of a variable. As shown in the statistical summary above, the variable with the highest coefficient of variation is volume, that represent the total sales value; instead, median_price has the lowest relative dispersion.

Regarding the asymmetry of the variables, this can be evaluated with Fisher’s skewness, a measure of deviation from symmetry. Once again, volume is the most asymmetric variable, as it exhibits the highest absolute value of skewness. Specifically, it shows positive skewness, meaning that most observations are concentrated on the left side of the distribution (at low values), while a few very high values extend the right tail of the distribution.

We can visualize all these characteristics of the volume distribution using a density plot.

# Density plot of volume
ggplot(dataset, aes(x = volume)) +
  geom_density(col = "black", fill = "skyblue") +
  labs(
    title = "Density plot of volume",
    x = "Volume (million $)",
    y = "Density"
  ) +
  scale_x_continuous(breaks = seq(0, max(volume), by = 10)) +
  coord_cartesian(xlim = c(min(volume), max(volume))) +
  theme_light() +
  theme(
    plot.title = element_text(size = 18, hjust = 0.5, face = "bold", color = "red3"),
    axis.text.x = element_text(face = "bold"),
    axis.text.y = element_text(face = "bold"),
    axis.title.x = element_text(size = 13, margin = margin(t = 10)),
    axis.title.y = element_text(size = 13, margin = margin(r = 10))
  )

4 Creating classes for a quantitative variable

The first goal of this section is to select the quantitative variable sales, divide it into classes and create a frequency distribution.
Let’s start computing the range of interval of sales.

# Check the min and max values of sales
min_sales <- stats["Minimum", "sales"]
max_sales <- stats["Maximum", "sales"]
range_sales <- stats["Range", "sales"]

cat("Sales minimum value:", min_sales, "\n")
cat("Sales maximum value:", max_sales, "\n")
cat("Sales range:", range_sales)
## Sales minimum value: 79 
## Sales maximum value: 423 
## Sales range: 344

We can see that sales has a minimum value of 79 and a maximum value of 423, with a range of 344. We may choose to divide the observations into \(7\) classes with amplitude \(50\), starting from \(75\) and ending at \(425\): in this way, we obtain classes that are easy to interpret and that cover the entire range of the sales variable. We now create these classes and compute the corresponding frequency distribution.

# Divide sales into 7 classes with amplitude 50
sales_classes <- cut(sales, breaks = seq(from=75, to=425, by=50))

# Evaluate frequencies
N <- length(sales)
ni_sales <- table(sales_classes)
fi_sales <- ni_sales/N
Ni_sales <- cumsum(ni_sales)
Fi_sales <- Ni_sales/N

# Create the dataframe of frequencies
sales_freq <- data.frame(
  labels = names(ni_sales),
  freq = as.vector(ni_sales),
  rel_freq = as.vector(fi_sales),
  cum_freq = as.vector(Ni_sales),
  cum_rel_freq = as.vector(Fi_sales)
)

colnames(sales_freq) <- c("Sales classes", "Frequency", "Relative frequency",
                          "Cumulative frequency", "Cumulative relative frequency")


# Display the frequency distribution rounded with 4 decimal
show_table(sales_freq, digits = 4)
Sales classes Frequency Relative frequency Cumulative frequency Cumulative relative frequency
(75,125] 59 0.2458 59 0.2458
(125,175] 61 0.2542 120 0.5000
(175,225] 46 0.1917 166 0.6917
(225,275] 30 0.1250 196 0.8167
(275,325] 26 0.1083 222 0.9250
(325,375] 13 0.0542 235 0.9792
(375,425] 5 0.0208 240 1.0000

We can observe that both the first quartile and the median value fall in the class \((125,175]\), while the third quartile falls in the class \((225,275]\). The modal class is \((125,175]\).

Now let’s visualize the computed frequencies using a column chart.

# Set the correct order of the classes for the plot

sales_freq$`Sales classes` <- factor(
  sales_freq$`Sales classes`,
  levels = sales_freq$`Sales classes`
)

# Plot a column chart of sales frequency
ggplot(sales_freq, aes(x = `Sales classes`, y = Frequency)) +
  geom_col(color = "black", fill = "darkorange") +
  geom_text(aes(label = Frequency), vjust = -0.3, size=3, fontface="bold") +
  labs(
    title = "Distribution of sales classes",
    x = "Sales classes",
    y = "Frequency"
  ) + 
  scale_y_continuous(breaks=seq(0, max(sales_freq$Frequency), by=10)) +
  theme_minimal() +
  coord_cartesian(ylim = c(0, max(sales_freq$Frequency) * 1.05)) +
  theme(
    plot.title = element_text(size=18, hjust = 0.5, face = "bold", color="red3"),
    axis.text.x = element_text(face = "bold"),
    axis.text.y = element_text(face = "bold"),
    axis.title.x = element_text(size = 13, margin = margin(t = 10)),
    axis.title.y = element_text(size = 13, margin = margin(r = 10))
    )

We can notice that the distribution of sales classes is right skewed, meaning that most records are concentrated at low sales class values.

Now we want to calculate the normalized Gini index to measure the heterogeneity of the frequency distribution across the classes of the sales variable. The normalized Gini index has a range of \([0, 1]\) and is computed as follows: \[ G = \frac{J}{J-1} \left(1 - \sum_{i=1}^{J} f_i^2 \right) \] where \(f_i\) are the relative frequencies of the \(J\) classes. Let’s create a function to calculate the normalized Gini index.

# ---------- NORMALIZED GINI INDEX FUNCTION ----------

gini_norm_index <- function(x){
  # x = variable to calculate normalized Gini index
  ni <- table(x)
  fi <- ni/length(x)
  J <- length(ni)
  gini_norm <- (J / (J-1)) * (1 - sum(fi^2))
  
  return(gini_norm)
}

# ----------------------------------------------------
# Evaluate normalized Gini index of sales classes
gini_norm_sales <- gini_norm_index(sales_classes)

cat("Normalized Gini index of sales classes =", round(gini_norm_sales, 4))
## Normalized Gini index of sales classes = 0.9421

The normalized Gini index is very close to \(1\), which means that the class frequencies of sales are highly heterogeneously distributed.

5 Probability computing

This section will present some examples of probability calculations, let’s look at them one by one. We will use the frequency distribution computed in the previous sections.

Let’s start by calculating the probability that a random observation in the dataset belongs to the city “Beaumont”: this corresponds to the relative frequency of that city in the cities frequency distribution.

# Probability that city is Beaumont
p_Beaumont <- freqs_city[freqs_city$city == "Beaumont", "Relative frequency"]
cat(paste("Probability (city = Beaumont) =", round(100*p_Beaumont, 2), "%"))
## Probability (city = Beaumont) = 25 %

Now let’s compute the probability that a random record contains the month of July: once again, this corresponds to the relative frequency of the month “Jul” in the monthly frequency distribution.

p_July <- freqs_month[freqs_month$month == "Jul", "Relative frequency"]
cat(paste("Probability (month = July) =", round(100*p_July, 2), "%"))
## Probability (month = July) = 8.33 %

Finally, let’s evaluate the probability that a random observation refers to the month of December 2012: it corresponds to the compound probability that a random record contains the month “Dec” and the year 2012, so it’s equal to the product of the two relative frequencies.

p_December <- freqs_month[freqs_month$month == "Dec", "Relative frequency"]
p_2012 <- freqs_year[freqs_year$year == "2012", "Relative frequency"]
p_December_2012 <- p_December * p_2012

cat(paste("Probability (date = December 2012) =", round(100*p_December_2012, 2), "%"))
## Probability (date = December 2012) = 1.67 %

6 Creating new variables

Let’s create, in the dataset, two new columns using the available variables:

  • A new column contains the average property price per unit sold.
    The mean price is computed as the volume (total sales value in million of dollars) divided by the sales (the total number of sales). Then, the mean_price column will be converted in dollars to align it with the unit of measurement for the median price.

  • A new column measures the efficiency of sales advertisements.
    A useful indicator for estimating the effectiveness of sales ads can be defined as the sales divided by the listings; this metric measures the average number of sales made per active listing and, therefore, how quickly the ads convert into sales.

# Evaluate the mean price, avoiding division by zero
dataset$mean_price <- ifelse(sales == 0, NA, (volume/sales) * 10^6)       

# Evaluate the efficiency of listings, avoiding division by zero
dataset$efficiency <- ifelse(listings == 0, NA, sales/listings)

attach(dataset)

# Show the dataset with the two new variables
show_table(head(dataset, 5))
city year month sales volume median_price listings months_inventory date mean_price efficiency
Beaumont 2010 Jan 83 14.16 163800 1533 9.5 2010-01-01 170626.5 0.05
Beaumont 2010 Feb 108 17.69 138200 1586 10.0 2010-02-01 163796.3 0.07
Beaumont 2010 Mar 182 28.70 122400 1689 10.6 2010-03-01 157697.8 0.11
Beaumont 2010 Apr 200 26.82 123200 1708 10.6 2010-04-01 134095.0 0.12
Beaumont 2010 May 202 28.83 123100 1771 10.9 2010-05-01 142737.6 0.11

We can summarize the statistical measure of the two new quantitative variables.

show_table(stat_summary(dataset, c("mean_price", "efficiency")))
mean_price efficiency
Mean value 154320.37 0.12
Minimum 97010.20 0.05
First quartile (Q1) 132938.94 0.09
Median value 156588.48 0.11
Third quartile (Q3) 173915.15 0.13
Maximum 213233.94 0.39
Variance 736984385.27 0.00
Standard deviation 27147.46 0.05
Range 116223.74 0.34
Interquartile range 40976.22 0.05
Coefficient of variation 0.18 0.39
Fisher’s Skewness -0.07 2.09
Excess Kurtosis -0.78 6.88

Let’s visualize the density plot of mean_price.

# Density plot of mean price
ggplot(dataset, aes(x = mean_price)) +
  geom_density(col = "black", fill = "skyblue") +
  labs(
    title = "Density plot of mean price",
    x = "Mean price ($)",
    y = "Density"
  ) +
  scale_x_continuous(breaks = seq(0, max(mean_price), by = 25000)) +
  coord_cartesian(xlim = c(min(mean_price), max(mean_price))) +
  theme_light() +
  theme(
    plot.title = element_text(size = 18, hjust = 0.5, face = "bold", color = "red3"),
    axis.text.x = element_text(face = "bold"),
    axis.text.y = element_text(face = "bold"),
    axis.title.x = element_text(size = 13, margin = margin(t = 10)),
    axis.title.y = element_text(size = 13, margin = margin(r = 10))
  )

Looking at the variable mean_price, we observe a a small coefficient of variation, indicating limited variability in the distribution; in addition, the distribution shows a slightly negative skewness and is platykurtic. We can assirt that mean_price appears to be a fairly stable indicator for measuring real estate prices.

Now let’s visualize the density plot of efficiency.

# Density plot of efficiency
ggplot(dataset, aes(x = efficiency)) +
  geom_density(col = "black", fill = "skyblue") +
  labs(
    title = "Density plot of efficiency",
    x = "Efficiency",
    y = "Density"
  ) +
  coord_cartesian(xlim = c(min(efficiency), max(efficiency))) +
  theme_light() +
  theme(
    plot.title = element_text(size = 18, hjust = 0.5, face = "bold", color = "red3"),
    axis.text.x = element_text(face = "bold"),
    axis.text.y = element_text(face = "bold"),
    axis.title.x = element_text(size = 13, margin = margin(t = 10)),
    axis.title.y = element_text(size = 13, margin = margin(r = 10))
  )

Regarding the variable efficiency, we observe a higher coefficient of variation, indicating moderate variability. Moreover, the distribution is strongly positively skewed and exhibits very high positive kurtosis, which means it is significantly leptokurtic.
The pronounced right tail tells us that most observations show low efficiency, while only a few observations are better at converting listings into sales.

Overall, the listing efficacy varies considerably across observations, suggesting that efficiency is not uniform throughout the dataset.

7 Conditional analysis

In this section we will perform conditional statistical analysis by city, year and month. More specifically, the goal is to consider the quantitative variables and calculate statistics for each category of the qualitative variables.

First, let’s create a function to perform a conditional analysis on a given set of quantitative variables, grouped by one categorical variable.

# -------------------- CONDITIONAL ANALYSIS FUNCTION --------------------

conditional_analysis <- function(df, quant_vars, categ_var){
  # df: data frame containing the data to analyze
  # quant_vars: vector of quantitative variables to summarise
  # categ_var: categorical variable used for grouping

  # Compute mean and standard deviation of quantitative variables 
  # grouped by each category of the qualitative variable
  stats <- df %>%
    group_by(across(all_of(categ_var))) %>%
    summarise(
      across(
        all_of(quant_vars),
        list(avg = mean, std = sd),
        .names = "{.fn}({.col})"
      ),
      .groups = "drop"
    )
  
  # Reshape the summary table for better visualization
  stats_transposed <- stats %>%
    pivot_longer(
      cols = -all_of(categ_var),
      names_to = "statistic",
      values_to = "value"
    ) %>%
    pivot_wider(
      names_from = all_of(categ_var),
      values_from = value
    )

  return(stats_transposed)
}

# -----------------------------------------------------------------------

We will apply this function to all the quantitative variables in the dataset, grouped by the categorical variables city, year, and month. We will report the mean and the standard deviation, and then plot bar charts showing the average volume, used here as the main variable of interest, for each categorical variable.

7.1 Analysis by city

Let’s summarize the statistics grouped by city, and then plot a bar chart of average volume by city.

# Summarize statistics grouping by city
city_stats <- conditional_analysis(
  dataset,
  c("sales", "volume", "median_price", "listings", "months_inventory", "mean_price", "efficiency"),
  "city"
)

# Show statistics rounded to 2 decimal places
show_table(city_stats)
statistic Beaumont Bryan-College Station Tyler Wichita Falls
avg(sales) 177.38 205.97 269.75 116.07
std(sales) 41.48 84.98 61.96 22.15
avg(volume) 26.13 38.19 45.77 13.93
std(volume) 6.97 17.25 13.11 3.24
avg(median_price) 129988.33 157488.33 141441.67 101743.33
std(median_price) 10104.99 8852.24 9336.54 11320.03
avg(listings) 1679.32 1458.13 2905.05 909.58
std(listings) 91.13 252.53 226.75 73.76
avg(months_inventory) 9.97 7.66 11.32 7.82
std(months_inventory) 1.65 2.25 1.89 0.78
avg(mean_price) 146640.41 183534.29 167676.76 119430.00
std(mean_price) 11232.13 15149.35 12350.51 11398.48
avg(efficiency) 0.11 0.15 0.09 0.13
std(efficiency) 0.03 0.07 0.02 0.02 
# Convert data frame in long format, filtering average volume
volume_by_city <- city_stats %>%
  filter(statistic == "avg(volume)") %>%
  pivot_longer(
    cols = -statistic,
    names_to = "city",
    values_to = "value"
  )

# Plot a bar chart of average volume by city
ggplot(volume_by_city, aes(
  x = reorder(city, value),
  y = value
  ))+
  geom_col(color = "black", fill = "darkseagreen") +
  labs(
    title = "Average volume by city",
    x = "City",
    y = "Average volume (million $)"
  ) + 
  scale_y_continuous(breaks=seq(0, max(volume_by_city$value), by=5)) +
  theme_minimal() +
  theme(
    plot.title = element_text(size = 18, hjust = 0.5, face = "bold", color = "red3"),
    axis.text.x = element_text(face = "bold"),
    axis.text.y = element_text(face = "bold"),
    axis.title.x = element_text(size = 13, margin = margin(t = 10)),
    axis.title.y = element_text(size = 12, margin = margin(r = 10))
  )

We can see that Tyler is the Texas city with the highest average volume, followed by Bryan-College Station and Beaumont, while Wichita Falls has the lowest average volume. From the summary above we can also notice that Bryan-College Station, followed by Tyler, shows very high dispersion in the volume distribution, while the variability of volume in Beaumont and Wichita Falls remains relatively low.

7.2 Analysis by year

Now we summarize the statistics grouped by year, and then plot a bar chart of average volume by year.

# Summarize statistics grouping by year
year_stats <- conditional_analysis(
  dataset,
  c("sales", "volume", "median_price", "listings", "months_inventory", "mean_price", "efficiency"),
  "year"
)

# Show statistics rounded to 2 decimal places
show_table(year_stats)
statistic 2010 2011 2012 2013 2014
avg(sales) 168.67 164.12 186.15 211.92 230.60
std(sales) 60.54 63.87 70.91 84.00 95.51
avg(volume) 25.68 25.16 29.27 35.15 39.77
std(volume) 10.80 12.20 14.52 17.93 21.19
avg(median_price) 130191.67 127854.17 130077.08 135722.92 139481.25
std(median_price) 21821.76 21317.80 21431.52 21708.08 25625.41
avg(listings) 1826.00 1849.65 1776.81 1677.60 1560.04
std(listings) 785.02 780.38 738.45 743.52 706.71
avg(months_inventory) 9.97 10.90 9.88 8.15 7.06
std(months_inventory) 2.08 2.07 1.61 1.69 1.75
avg(mean_price) 150188.58 148250.63 150898.68 158705.25 163558.70
std(mean_price) 23279.55 24938.38 26438.50 26523.81 31740.53
avg(efficiency) 0.10 0.09 0.11 0.13 0.16
std(efficiency) 0.03 0.02 0.03 0.04 0.06 
# Convert data frame in long format, filtering average volume
volume_by_year <- year_stats %>%
  filter(statistic == "avg(volume)") %>%
  pivot_longer(
    cols = -statistic,
    names_to = "year",
    values_to = "value"
  )

# Plot a bar chart of average volume by year
ggplot(volume_by_year, aes(
  x = year,
  y = value
  ))+
  geom_col(color = "black", fill = "wheat") +
  labs(
    title = "Average volume by year",
    x = "Year",
    y = "Average volume (million $)"
  ) + 
  scale_y_continuous(breaks=seq(0, max(volume_by_year$value), by=5)) +
  theme_minimal() +
  theme(
    plot.title = element_text(size = 18, hjust = 0.5, face = "bold", color = "red3"),
    axis.text.x = element_text(face = "bold"),
    axis.text.y = element_text(face = "bold"),
    axis.title.x = element_text(size = 13, margin = margin(t = 10)),
    axis.title.y = element_text(size = 12, margin = margin(r = 10))
  )

From the bar chart above we can see that the total average volume of real estate was approximately \(26\) million dollars in 2010, then it decreased slightly in 2011, and gradually increased in the following years, reaching about \(40\) million dollars in 2014. We can also observe that this increase is not driven by significant changes in the mean property price, but rather by an increase in sales over the years. Moreover, we can observe that volume variability has generally increased over the years, more than doubling between 2010 and 2014.

7.3 Analysis by month

Finally, we summarize the statistics grouped by month, and then plot a bar chart of average volume by month.

# Summarize statistics grouping by month
month_stats <- conditional_analysis(
  dataset,
  c("sales", "volume", "median_price", "listings", "months_inventory", "mean_price", "efficiency"),
  "month"
)

# Show statistics rounded to 2 decimal places
show_table(month_stats, scrollable = TRUE)
statistic Jan Feb Mar Apr May Jun Jul Aug Sep Oct Nov Dec
avg(sales) 127.40 140.85 189.45 211.70 238.85 243.55 235.75 231.45 182.35 179.90 156.85 169.40
std(sales) 43.38 51.07 59.18 65.40 83.12 95.00 96.27 79.23 72.52 74.95 55.47 60.75
avg(volume) 19.00 21.65 29.38 33.30 39.70 41.30 39.12 38.01 29.60 29.08 24.81 27.09
std(volume) 8.37 10.09 12.02 14.52 19.02 21.08 21.41 18.05 15.22 15.13 11.15 12.57
avg(median_price) 124250.00 130075.00 127415.00 131490.00 134485.00 137620.00 134750.00 136675.00 134040.00 133480.00 134305.00 133400.00
std(median_price) 25151.28 22822.59 23442.03 21458.40 18796.26 19231.02 21944.78 22488.38 24344.10 26358.07 24691.47 22809.76
avg(listings) 1647.05 1692.50 1756.70 1825.70 1823.85 1833.25 1821.20 1786.30 1748.90 1710.35 1652.70 1557.75
std(listings) 704.61 711.20 727.35 770.43 790.22 811.63 826.72 815.87 802.66 779.16 741.25 692.57
avg(months_inventory) 8.84 9.06 9.40 9.72 9.68 9.70 9.62 9.39 9.19 8.94 8.66 8.12
std(months_inventory) 1.97 1.98 2.06 2.24 2.38 2.41 2.50 2.45 2.52 2.44 2.37 2.27
avg(mean_price) 145640.42 148840.48 151136.54 151461.33 158235.03 161545.82 156881.00 156455.56 156522.32 155897.37 154233.00 154995.52
std(mean_price) 29819.11 25120.42 23237.92 26174.30 25787.19 23470.46 27220.12 28253.21 29669.41 32527.29 29684.87 27008.87
avg(efficiency) 0.08 0.09 0.12 0.13 0.14 0.14 0.14 0.14 0.11 0.11 0.10 0.12
std(efficiency) 0.02 0.02 0.03 0.04 0.05 0.06 0.07 0.05 0.03 0.04 0.03 0.04 
# Convert data frame in long format, filtering average volume
volume_by_month <- month_stats %>%
  filter(statistic == "avg(volume)") %>%
  pivot_longer(
    cols = -statistic,
    names_to = "month",
    values_to = "value"
  )

# Order months in calendar order
volume_by_month$month <- factor(
  volume_by_month$month,
  levels = c("Jan", "Feb", "Mar", "Apr", "May", "Jun",
             "Jul", "Aug", "Sep", "Oct", "Nov", "Dec")
)

# Plot a bar chart of sales by month
ggplot(volume_by_month, aes(
  x = month,
  y = value
  ))+
  geom_col(color = "black", fill = "lightsteelblue") +
  labs(
    title = "Average volume by month",
    x = "Month",
    y = "Average volume (million $)"
  ) + 
  scale_y_continuous(breaks=seq(0, max(volume_by_month$value), by=5)) +
  theme_minimal() +
  theme(
    plot.title = element_text(size = 18, hjust = 0.5, face = "bold", color = "red3"),
    axis.text.x = element_text(face = "bold"),
    axis.text.y = element_text(face = "bold"),
    axis.title.x = element_text(size = 13, margin = margin(t = 10)),
    axis.title.y = element_text(size = 12, margin = margin(r = 10))
  )

This chart is useful for visualizing changes in average volume across the different months, regardless of the reference year. It allows us to assess the seasonality of volume: we can observe that June has the highest average volume, followed by May and July; on the other hand, January has the lowest average volume, followed by February and November. Moreover, the period from May to August shows the highest average volume, due to an increase in sales. We can also observe that volume dispersion generally follows the same trend of average monthly volume, so the months with the highest volume are also those with the highest volume variability, and vice versa.

8 Creating visualizations with ggplot2

To enhance our statistical analysis and gain a deeper understanding of the real estate market, this section presents a variety of visualizations created with the ggplot2 package.

8.1 Boxplot of median price by city

Let’s visualize a boxplot showing the distribution of median_price by city.

# Boxplot of median price by city
ggplot(dataset, aes(
  x = city,
  y = median_price
  )) +
  geom_boxplot(
    color = "black",
    fill = "mistyrose",
    outlier.size = 2
  ) +
  labs(
    title = "Boxplot of median price by city",
    x = "City",
    y = "Median price ($)"
  ) + 
  theme_light() +
  theme(
    plot.title = element_text(size = 18, hjust = 0.5, face = "bold", color = "red3"),
    axis.text.x = element_text(face = "bold"),
    axis.text.y = element_text(face = "bold"),
    axis.title.x = element_text(size = 13, margin = margin(t = 10)),
    axis.title.y = element_text(size = 13, margin = margin(r = 10))
  )

We can summarize the boxplot statistics in the table below, including the number of outliers for each city, and the skewness of median price distribution.

# Boxplot statistics of median price by city
outlier_median_price_by_city <- dataset %>%
  group_by(city) %>%
  summarise(
    Q1 = quantile(median_price, 0.25),
    median = median(median_price),
    Q3 = quantile(median_price, 0.75),
    IQR = IQR(median_price),
    lower_bound = Q1 - 1.5 * IQR,
    upper_bound = Q3 + 1.5 * IQR,
    n_outliers = sum(median_price < lower_bound | median_price > upper_bound),
    skewness = skewness(median_price)
  )

show_table(outlier_median_price_by_city)
city Q1 median Q3 IQR lower_bound upper_bound n_outliers skewness
Beaumont 123025 130750 134550 11525 105737.5 151837.5 1 0.36
Bryan-College Station 150800 155400 161975 11175 134037.5 178737.5 1 0.71
Tyler 133975 142200 147675 13700 113425.0 168225.0 0 0.12
Wichita Falls 92800 102300 109175 16375 68237.5 133737.5 1 0.22

We can extract some information from the boxplot and the summary:

  • Bryan-College Station has the highest median property price and the smallest interquartile range, indicating the lowest variability in median price.
  • Wichita Falls has the lowest real estate median price and shows the greatest dispersion.
  • Beaumont, Bryan-College Station and Wichita Falls each have one outlier bove the upper bound (\(Q_3 + 1.5 \cdot IQR\)), while Tyler is the only city without outliers.
  • All cities have a slighly right-skewed distribution of median price. This suggests that most property median prices are concentrated at lower values, while a few higher values stretch the distribution to the right.

8.2 Bar charts of sales by month and city

Now we present a bar chart illustrating the total sales by month and city. This visualization will help us to assess the seasonality of sales across the different cities.

# Compute total sales by city and month
sales_by_city_month <- dataset %>%
  group_by(city, month) %>%
  summarise(total_sales = sum(sales))

# Plot a bar chart of sales by city and month
ggplot(sales_by_city_month, aes(
  x = month,
  y = total_sales,
  fill = city
  ))+
  geom_col(
    position = position_dodge(width = 0.9),
    width = 0.6
    ) +
  labs(
    title = "Total sales by month and city",
    x = "Month",
    y = "Sales",
    fill = "City"
  ) + 
  theme_minimal() +
  theme(
    plot.title = element_text(size = 18, hjust = 0.5, face = "bold", color = "red3"),
    axis.text.x = element_text(face = "bold"),
    axis.text.y = element_text(face = "bold"),
    axis.title.x = element_text(size = 13, margin = margin(t = 10)),
    axis.title.y = element_text(size = 13, margin = margin(r = 10))
  )

From the bar chart above we can observe that:

  • Comparing the different cities, Tyler has the highest number of sales across all months, while Wichita Falls has the lowest.
  • All cities have their lowest total sales in January.
  • Tyler and Bryan-College Station have their highest total sales in June, while Beaumont reaches its maximum sales in August and Wichita Falls in May.

Now let’s plot a normalized stacked bar chart showing the proportion of sales by month and city.

# Plot a normalized stacked bar chart of sales by city and month
ggplot(sales_by_city_month, aes(
  x = month,
  y = total_sales,
  fill = city
  ))+
  geom_col(position = "fill", width = 0.6) +
  scale_y_continuous(
    labels = scales::percent,
    breaks = seq(0, 1, 0.1)
  ) +
  labs(
    title = "Percentage of sales by month and city",
    x = "Month",
    y = "Percentage of Sales",
    fill = "City"
  ) + 
  theme_minimal() +
  theme(
    plot.title = element_text(size = 18, hjust = 0.5, face = "bold", color = "red3"),
    axis.text.x = element_text(face = "bold"),
    axis.text.y = element_text(face = "bold"),
    axis.title.x = element_text(size = 13, margin = margin(t = 10)),
    axis.title.y = element_text(size = 12, margin = margin(r = 10))
  )

Once again, we can see that Tyler consistently represents the largest share of sales, followed by Bryan-College Station and Beaumont, while Wichita Falls accounts for the smallest percentage of sales. The distribution among the cities appears fairly stable over the months, with small fluctuations in each city’s share.

To the previous visualization, we can add the variable year showing one stacked bar chart for each year in the dataset.

# Plot normalized stacked bar charts of sales by city and month, one bar chart per year
ggplot(dataset, aes(
  x = month,
  y = sales,
  fill = city
  ))+
  geom_col(position = "fill", width = 0.6) +
  facet_wrap(~ year, ncol = 1) + 
  scale_x_discrete(drop = FALSE) +
  scale_y_continuous(
    labels = scales::percent,
    breaks = seq(0, 1, 0.2)
  ) +
  labs(
    title = "Percentage of sales by year, month and city",
    x = "Month",
    y = "Percentage of Sales",
    fill = "City"
  ) + 
  theme_minimal() +
  theme(
    plot.title = element_text(size = 22, hjust = 0.5, face = "bold", color = "red3"),
    strip.text = element_text(size = 14, face = "bold"),
    strip.background = element_rect(fill = "grey90", color = "grey50"),
    axis.text.x = element_text(size = 12, face = "bold"),
    axis.text.y = element_text(size = 10, face = "bold"),
    axis.title.x = element_text(size = 16, margin = margin(t = 10)),
    axis.title.y = element_text(size = 16, margin = margin(r = 10)),
    legend.title = element_text(size = 14, face = "bold"),
    legend.text = element_text(size = 12),
    legend.key.size = unit(0.8, "cm")
  )

The percentage composition of sales by month and city appears relatively stable across years. Tyler systematically accounts for the largest share of sales over years, while Wichita Falls has the lowest percentage of sales over years.

8.3 Line chart of sales trend

Let’s plot a line chart showing the sales trend using the date variable. This will allow us to observe how sales evolve across different time periods.

# Compute monthly sales
monthly_sales <- dataset %>%
  group_by(date) %>%
  summarise(sum_sales = sum(sales), .groups = "drop") %>%
  arrange(date)

# Line chart of sales trend
ggplot(monthly_sales, aes(
  x = date,
  y = sum_sales
  )) +
  geom_line(col = "lightblue3", linewidth = 1) +
  geom_point(col="steelblue4", size = 1.5) +
  scale_x_date(
    date_breaks = "3 months",
    date_labels = "%b %Y"
  ) +
  labs(
    title = "Sales trend",
    x = "Date",
    y = "Total Sales"
  ) + 
  theme_minimal() +
  theme(
    plot.title = element_text(size = 18, hjust = 0.5, face = "bold", color = "red3"),
    axis.text.x = element_text(size = 6, face = "bold"),
    axis.text.y = element_text(face = "bold"),
    axis.title.x = element_text(size = 13, margin = margin(t = 10)),
    axis.title.y = element_text(size = 13, margin = margin(r = 10))
  )

From the line chart above we can clearly observe a recurrent pattern in sales trend. Similar to the average monthly volume studied before, sales also exhibit a seasonal pattern, reaching their highest values in summer and their lowest in winter. Moreover, as with the average volume, the annual mean of sales decreases from 2010 to 2011 and then increases from 2011 to 2014.

8.4 Boxplot of volume by city and year

Let’s use a boxplot to visualize the distribution of volume grouped by city and year.

# Boxplot of volume by city and year
ggplot(dataset, aes(
  x = city,
  y = volume,
  fill = factor(year)
  )) +
  geom_boxplot(
    color = "black",
    outlier.size = 2
  ) +
  scale_fill_brewer(palette = "Set2") +
  labs(
    title = "Boxplot of volume by city and year",
    x = "City",
    y = "Volume (million $)",
    fill = "Year"
  ) + 
  theme_light() +
  theme(
    plot.title = element_text(size = 18, hjust = 0.5, face = "bold", color = "red3"),
    axis.text.x = element_text(face = "bold"),
    axis.text.y = element_text(face = "bold"),
    axis.title.x = element_text(size = 13, margin = margin(t = 10)),
    axis.title.y = element_text(size = 12, margin = margin(r = 10))
  )

From the boxplot we can extract some insights:

  • Wichita Falls is the city with the lowest median volume over the years, and it also shows the smallest variability.
  • Beaumont has a higher median volume than Wichita Falls, although values are still relatively low and the dispersion is limited. In 2012, it shows the lowest variability and three outliers.
  • Bryan-College Station has a higher median volume and the largest dispersion overall.
  • Tyler has the highest median volume across the years, with variability lower than that of Bryan-College Station but still high.

Overall, the boxplot suggests a growth trend in volume over time for Texas cities, with increasing dispersion in the cities with the highest volumes as seen before.

8.5 Line chart of efficiency by city and year

Now we want to plot a line chart showing the average efficiency (the mean of ratio sales/listings) grouped by city and year.
This will allow us to study how the advertisement efficiency varies over time and across cities.

# Compute the average efficiency by city and year
efficiency_by_city_year <- dataset %>%
  group_by(city, year) %>%
  summarise(avg_efficiency = mean(efficiency))

# Line chart of average efficiency by city and year
ggplot(efficiency_by_city_year, aes(
  x = year,
  y = avg_efficiency,
  color = city,
  group = city
  )) +
  geom_line(linewidth = 1.2) +
  geom_point(size = 1.8) +
  scale_y_continuous(breaks = seq(0, max(efficiency_by_city_year$avg_efficiency)*1.1, by = 0.025)) + 
  labs(
    title = "Average efficiency by city and year",
    x = "Year",
    y = "Average efficiency",
    color = "City"
  ) + 
  theme_minimal() +
  theme(
    plot.title = element_text(size = 18, hjust = 0.5, face = "bold", color = "red3"),
    axis.text.x = element_text(size = 6, face = "bold"),
    axis.text.y = element_text(face = "bold"),
    axis.title.x = element_text(size = 13, margin = margin(t = 10)),
    axis.title.y = element_text(size = 13, margin = margin(r = 10))
  )

From the line chart of average efficiency we can observe that:

  • Tyler is the only city with a monotonic increasing trend over time, although it remains the least efficient among cities.
  • Beaumont and Bryan-College Station both experienced a slight decline of average efficiency from 2010 to 2011, followed by an increase in the next years. The growth was modest in Beaumont, while Bryan-College Station showed a stronger upward trend, which made it the city with the highest efficiency overall in 2013 and 2014.
  • Wichita Falls, instead, followed a more unstable trend, with a decline of average efficiency from 2010 to 2011, a slight increase until 2013, and another decline in 2014.

9 Conclusions

From the statistical analysis conducted on the dataset from Texas Realty Insights, we can extract some useful information to understand the real estate market situation and optimize sales strategies. Let’s summarize the key insights.

  • General upward trend of market over the years.

    Annual sales and volume declined between 2010 and 2011, and then increase from 2011 to 2014.
    Advertisements efficiency followed a similar pattern,suggesting a potential association between improved listing effectiveness and the overall growth of the market over time.
    More specifically, total sales value was approximately 25.7 million $ in 2010 and grew to 39.8 million $ in 2014. The ads efficiency increased from 9.97% in 2010 to 15.7% in 2014, a good change that still has plenty of room for improvement.

  • Significant Seasonality of sales and volume.

    Both sales and volume exhibit a clear seasonal pattern, reaching their highest values during the summer and their lowest values in winter. June achieves the best performance in terms of sales and volume, whereas January shows the worst results. The company could take advantage of this seasonal trend by improving listing effectiveness during the winter months and maximizing sales during the summer.

  • Sales and volume vary considerably across cities.

    The ranking of cities in terms of both sales and volume is as follows:

    1. Tyler
    2. Bryan-College Station
    3. Beaumont
    4. Wichita Falls

    It has been observed that the percentage share of sales across cities has remained relatively stable over time.

  • Geographic variations in median prices.

    Median prices exhibit marked differences across cities. In descending order, the median price rank is as follows:

    1. Bryan-College Station
    2. Tyler
    3. Beaumont
    4. Wichita Falls

    It may be useful for the company to examine if there is a correlation between median property prices and the performance of each city.

  • Different trends in advertisements efficiency across cities.

    Tyler is the only city that shows a steady increase in ads efficiency over time. However, despite this upward trend, this metric remains the lowest among the cities. Given that Tyler also has the highest volume, improving advertisement efficiency could have a substantial impact on market performance.

    Beaumont and Bryan-College Station both experienced a slight decline in efficiency from 2010 to 2011, which could help to explain the modest decrease in sales and volume observed during that period. In the following years, from 2011 to 2014, listing efficacy increased, especially in Bryan-College Station, which is the city with the best advertisement efficiency overall.

    Wichita Falls, by contrast, exhibits a less stable trend, with a decline in efficiency from 2010 to 2011, a slight increase until 2013, and another decline in 2014. It may be interesting for the company to investigate the reasons behind this pattern.