在 R 中绘制双变量分布热力图

作为一名数据科学家，您需要分析数据集中各特征的分布。通常我们会使用直方图，它非常有助于展示变量的取值范围、偏离程度以及数值的集中区间。

本教程使用 R。如果您尚不熟悉 R，建议先在 DataCamp 上学习我们的免费课程Introduction to R。

共享单车数据集

在本教程中，您将使用来自 Kaggle 竞赛的Bike Sharing Demand Dataset。

打开链接进入数据页签，然后下载 train.csv 文件。

将数据集加载到 R 中。

library(dplyr)
bike<- read.csv("train.csv",
                na.strings = FALSE,
                strip.white = TRUE)
glimpse(bike)技术支持

数据集特征

特征选择

您将按 season 分组，分析 atemp 与 humidity 之间的关系，因此从数据集中选取这几个变量。记得加载 dplyr 包。

bike_data<-
  bike %>%  
  dplyr::select(season, atemp, humidity)
head(bike_data)技术支持

直方图

现在可以开始分析数据了，先为每个季节的各个特征创建直方图。我们将使用 ggplot2 库，请别忘了加载它。

湿度直方图

library(ggplot2)
bike_data %>%
  ggplot( aes(x=humidity) ) +
  geom_histogram(bins=30) +
  facet_wrap(~season,ncol = 2)技术支持

从直方图可以看到，冬季湿度更高（这很直观），同时也能感受到各季节湿度的分布形态。

现在加入密度曲线、密度地毯图，并修改分面名称，让可视化更美观。

bike_data %>%
  mutate(season_label = case_when(
    season == 1 ~ "Spring",
    season == 2 ~ "Summer",
    season == 3 ~ "fall",
    season == 4 ~ "winter")) %>%
  ggplot( aes(x=humidity) ) +
  geom_histogram(bins=30,aes(y = ..density..)) +
  geom_rug()+
  geom_density()+
  facet_wrap(~season_label,ncol = 2)+
  ggtitle("Histogram of humidity by season")技术支持

现在，您能更直观地把握湿度的变化规律了，用同样的方法处理 atemp。

atemp 直方图

复用上述代码，把变量改为 atemp：

bike_data %>%
  mutate(season_label = case_when(
    season == 1 ~ "Spring",
    season == 2 ~ "Summer",
    season == 3 ~ "fall",
    season == 4 ~ "winter")) %>%
  ggplot( aes(x=atemp) ) +
  geom_histogram(bins=30,aes(y = ..density..)) +
  geom_rug()+
  geom_density()+
  facet_wrap(~season_label,ncol = 2)+
  ggtitle("Histogram of atemp by season")技术支持

按季节划分的 atemp 直方图

目前为止一切顺利，各季节的图形都与预期一致。每张图都给出了概率密度函数的近似。

但仅凭这些图还很难直观地看到特征之间如何相互作用。首先遇到的困难是，您需要一个三维直方图，因为要寻找两个特征的联合分布函数。

为数据拟合双变量分布

目标是可视化双变量分布。为此，首先需要对数据拟合一个双变量分布。这里将使用 MASS 库和 kde2d 函数。kde2d 函数会在假设随机变量服从正态性的前提下估计双变量分布。该函数的输入如下：

现在将其应用于夏季数据：

library(MASS)
bike_data_summer <- bike_data %>% filter(season==1)
bike_density <- kde2d(bike_data_summer$atemp,bike_data_summer$humidity, n=1000)技术支持

注意： 加载 MASS 库时会带有一个 select 函数，这会与 dplyr 的 select 产生冲突。为解决该问题，使用 dplyr 的 select 时需加上命名空间前缀 dplyr::，例如 dplyr::select()，以确保正常工作。

现在查看 bike_density 的类与结构：

class(bike_density)
str(bike_density)技术支持

它是一个包含 x、y 和 z 三个值的列表。前两个是特征取值，第三个是一个矩阵，包含概率密度函数（pdf）的数值。 使用 contour 函数查看密度：

contour(bike_density)
text(12.2,92,"1",cex=1.5)
points(12.2,89,col="red",pch=18)
text(10.4,47,"2",cex=1.5)
points(10.4,43,col="red",pch=18)
text(21,45,"3",cex=1.5)
points(21,40,col="red",pch=18)
text(21.6,82,"4",cex=1.5)
points(21.6,78,col="red",pch=18)技术支持

从图中可以识别出一些累积点，可理解为这两个特征组合中最常出现的取值。各点的数值汇总如下表：

热力图

要创建热力图，首先需要定义要使用的颜色刻度，可通过 colorRampPalette 函数完成。

hm_col_scale<-colorRampPalette(c("black","blue","green","orange","red"))(1000)
技术支持

hm_col_scale 是一个颜色向量，从黑到红，以 RGB 颜色在 1000 个步长中渐变。

现在使用 image 函数绘制热力图，该函数用于绘制矩阵。

image(bike_density$z,  
      col = hm_col_scale,
      zlim=c(min(bike_density$z), max(bike_density$z)))
text(0.31,0.46,"Most Frequent",col="blue",cex=0.8)
points(0.31,0.45,col="blue",pch=18)
技术支持

如您所见，相较于等高线图，这里包含了更多信息，因为颜色展示了累积点的强度。可以看到，春季湿度与 atemp 最常见的取值大致在点 $(10.4,43)$ 附近。

总结

现在您已经学会了如何创建热力图，接下来就轮到您自己动手了。先去看看其他季节以及它们与其他变量之间的关系。

思考如何利用这些信息，在不同变量取值下，根据最常见点具有更高概率这一事实，来判断需要投放多少辆单车。

如果您想进一步学习如何使用 ggplot 生成数据可视化，欢迎查看我们的课程 Data Visualization with ggplot2。