探索性數據分析（EDA）
探索性數據分析exploratory data analysis
1 對分布進行可視化表示
分類變量在 R 中通常保存為因子或字符向量。要想檢查分類變量的分布，可以使用條形圖：
ggplot(data = diamonds) + geom_bar(mapping = aes(x = cut))
條形的高度表示每個 x 值中觀測的數量，你可以使用 dplyr::count() 手動計算出這些值：
diamonds %>% count(cut)
要想檢查連續變量的分布，可以使用直方圖：
ggplot(data = diamonds) + geom_histogram(mapping = aes(x = carat), binwidth = 0.5)
最高的條形表示幾乎有 30 000 個觀測的 carat 值在 0.25 和 0.75 之間，這兩個值分別是條形的左側值和右側值。可以通過 dplyr::count() 和 ggplot2::cut_width() 函數的組合來手動計算結果：
diamonds %>% count(cut_width(carat, 0.5))
如果只考慮重量小於 3 克拉的鑽石，並選擇一個更小的分箱寬度，那么直方圖如下所示：
smaller <- diamonds %>% filter(carat < 3) ggplot(data = smaller, mapping = aes(x = carat)) + geom_histogram(binwidth = 0.1)
如果想要在同一張圖上疊加多個直方圖，那么我們建議你使用 geom_freqploy() 函數來代替
geom_histogram() 函數。 ggplot(data = smaller, mapping = aes(x = carat, color = cut)) + geom_freqpoly(binwidth = 0.1)
1.2 典型值
以下的直方圖顯示了美國黃石國家公園中的老忠實噴泉的 272 次噴發的時長（單位為分
鍾）。噴發時間似乎聚集成了兩組：短噴發（2 分鍾左右）和長噴發（4~5 分鍾），這兩組
間幾乎沒有其他噴發時間：
ggplot(data = faithful, mapping = aes(x = eruptions)) + geom_histogram(binwidth = 0.25)
1.3 異常值
查看鑽石數據集中 y 軸變量的分布，唯一能表示存在異常值的證據是，y 軸的取值范圍出奇得寬：
ggplot(diamonds) + geom_histogram(mapping = aes(x = y), binwidth = 0.5)
正常值分箱中的觀測太多了，以致於包括異常值的分箱高度太低，因此我們根本看不見
（如果仔細觀察 x 軸 0 刻度附近，沒准你能發現點什么）。為了更容易發現異常值，我們可以使用 coord_cartesian() 函數將 y 軸靠近 0 的部分放大：
ggplot(diamonds) + geom_histogram(mapping = aes(x = y), binwidth = 0.5) + coord_cartesian(ylim = c(0, 50))
（coord_cartesian() 函數中有一個用於放大 x 軸的 xlim() 參數。我們就可以看出有 3 個異常值，分別位於 0、30 左右和 60 左右。我們使用dplyr 將它們找出來：
unusual <- diamonds %>% filter(y < 3 | y > 20) %>% arrange(y) unusual
結果分析：y 變量測量鑽石的三個維度之一，單位為毫米。我們知道鑽石的寬度不可能是 0 毫米，因此這些值肯定是錯誤的。我們也完全可以認為 32 毫米和 59 毫米同樣是令人難以置信的，這樣的鑽石長度超過 1 英寸（1 英寸 =2.54 厘米），簡直就是無價之寶！
1.4 缺失值
如果在數據集中發現異常值，將帶有可疑值的行全部丟棄：
diamonds2 <- diamonds %>% filter(between(y, 3, 20))
使用缺失值來代替異常值。最簡單的做法就是使用 mutate() 函數創建一個新變量來代替原來的變量。你可以使用 ifelse() 函數將異常值替換為 NA：
diamonds2 <- diamonds %>% mutate(y = ifelse(y < 3 | y > 20, NA, y))
ggplot2 也遵循不能無視缺失值的原則。所以ggplot2 在繪圖時會忽略缺失值，但會提出警告以通知缺失值被丟棄了：
ggplot(data = diamonds2, mapping = aes(x = x, y = y)) + geom_point()
要想不顯示這條警告，可以設置 na.rm = TRUE：
ggplot(data = diamonds2, mapping = aes(x = x, y = y)) + geom_point(na.rm = TRUE)
在nycflights13::flights 中，dep_time 變量中的缺失值表示航班取消了。因此，你應該比較
一下已取消航班和未取消航班的計划出發時間。可以使用 is.na() 函數創建一個新變量來完成這個操作：
nycflights13::flights %>%
mutate(
cancelled = is.na(dep_time),
sched_hour = sched_dep_time %/% 100,
sched_min = sched_dep_time %% 100,
sched_dep_time = sched_hour + sched_min / 60
) %>%
ggplot(mapping = aes(sched_dep_time)) +
geom_freqpoly(
mapping = aes(color = cancelled),
binwidth = 1/4
)
2 相關變動
2.1 分組變量與連續變量
我們探索一下鑽石價格是如何隨着質量而變化的：
ggplot(data = diamonds, mapping = aes(x = price)) + geom_freqpoly(mapping = aes(color = cut), binwidth = 500)
為了讓比較變得更容易，需要改變 y 軸的顯示內容，不再顯示計數，而是顯示密度。密度
是對計數的標准化，這樣每個頻率多邊形下邊的面積都是 1：
ggplot(
data = diamonds,
mapping = aes(x = price, y = ..density..)
) +
geom_freqpoly(mapping = aes(color = cut), binwidth = 500)
結果分析：這張圖的部分內容非常令人驚訝，其顯示出一般鑽石（質量最差）的平均價格是最高的！但這可能是因為頻率多邊形圖很難解釋，所以這張圖還有很多可以改進的地方。
箱線圖是對變量值分布的一種簡單可視化表示，這種圖在統計學家中非常流行。每張箱線圖都包括以下內容。
使用 geom_boxplot() 函數查看按切割質量分類的價格分布：
ggplot(data = diamonds, mapping = aes(x = cut, y = price)) + geom_boxplot()
cut 是一個有序因子：“一般”不如“較好”、“較好”不如“很好”，以此類推。因為很多分類變量並沒有這種內在的順序，所以有時需要對其重新排序來繪制信息更豐富的圖形。重新排序的其中一種方法是使用 reorder() 函數。
我們看一下 mpg 數據集中的 class 變量。你可能很想知道公路里程因汽車類別的不同會有怎樣的變化：
ggplot(data = mpg, mapping = aes(x = class, y = hwy)) + geom_boxplot()
基於 hwy 值的中位數對 class 進行重新排序：
ggplot(data = mpg) +
geom_boxplot(
mapping = aes(
x = reorder(class, hwy, FUN = median),
y = hwy
)
)
如果變量名很長，那么將圖形旋轉 90 度效果會更好一些。你可以通過 coord_flip() 函數
完成這一操作：
ggplot(data = mpg) +
geom_boxplot(
mapping = aes(
x = reorder(class, hwy, FUN = median),
y = hwy
)
) +
coord_flip()
2.2 兩個分類變量
要想對兩個分類變量間的相關變動進行可視化表示，需要計算出每個變量組合中的觀測數
量。完成這個任務的其中一種方法是使用內置的 geom_count() 函數：
ggplot(data = diamonds) + geom_count(mapping = aes(x = cut, y = color))
圖中每個圓點的大小表示每個變量組合中的觀測數量。
計算變量組合中的觀測數量的另一種方法是使用 dplyr：
diamonds %>% count(color, cut)
接着使用 geom_tile() 函數和填充圖形屬性進行可視化表示：
diamonds %>%
count(color, cut) %>%
ggplot(mapping = aes(x = color, y = cut)) +
geom_tile(mapping = aes(fill = n))
2.3 兩個連續變量
可以將相關變動看作點的模式。你可以看到鑽石的克拉數和價值之間存在一種指數關系：
ggplot(data = diamonds) + geom_point(mapping = aes(x = carat, y = price), alpha = 1 / 100) #使用 alpha 圖形屬性添加透明度
geom_bin2d() 和 geom_hex() 函數將坐標平面分為二維分箱，並使用一種填充顏色表示落入
每個分箱的數據點。geom_bin2d() 創建長方形分箱。geom_hex() 創建六邊形分箱。要想使
用 geom_hex()，需要安裝 hexbin 包：
ggplot(data = smaller) + geom_bin2d(mapping = aes(x = carat, y = price))
library(hexbin) ggplot(data = smaller) + geom_hex(mapping = aes(x = carat, y = price))
以對 carat 進行分箱，然后為每個組生成一個箱線圖：
ggplot(data = smaller, mapping = aes(x = carat, y = price)) + geom_boxplot(mapping = aes(group = cut_width(carat, 0.1)))
3 模式和模型
的美國黃石國家公園中老忠實噴泉的噴發時長和兩次噴發之間的等待時間做出一張散點圖，該圖會顯示出一個模式：較長的等待時間與較長的噴發時間是相關的。圖中還顯示出兩個簇：
ggplot(data = faithful) + geom_point(mapping = aes(x = eruptions, y = waiting))
以下代碼擬合了一個模型，可以根據 carat 預測 price，並計算出殘差（預測值和實際值之間的差別）。一旦去除克拉數對價格的影響，殘差就能反映出鑽石的價格：
library(modelr) mod <- lm(log(price) ~ log(carat), data = diamonds) diamonds2 <- diamonds %>% add_residuals(mod) %>% mutate(resid = exp(resid)) ggplot(data = diamonds2) + geom_point(mapping = aes(x = carat, y = resid))
