R语言学习：图例管理，SQL语句操作数据框，tidyverse包实操，情感分析，词云图

library(dplyr)
library(nycflights13)
library(knitr)

# 探索数据框的方法
View(flights)
glimpse(flights)
kable(airlines)

library(ggplot2)
ToothGrowth$dose<-factor(ToothGrowth$dose)

# 图例保持
# 1）图的外面
p <- ggplot(ToothGrowth, aes(x=dose, y=len, fill=dose)) +
  geom_boxplot()
p
p +theme(legend.position = "bottom")
# 2）图的里面
p + theme(legend.position = c(.9, .9))
p + theme(legend.position = c(.9, .1))

# 图例移除
p + theme(legend.position = "none")

library(sqldf)

# R语言自带数据集mtcars
# 第一步：根据数据处理逻辑编辑SQL语句
sql_st = "
select 
*
from mtcars
where carb = 1
order by mpg 
"
# 第二步：使用sqldf包的sqldf函数执行SQL语句
new_df <- sqldf(sql_st, row.names = TRUE)

# 第三步：结果查看
new_df

library(tidyverse)
library(lubridate)
library(nycflights13)

# 获取数据
# 航班数据
head(flights)

# 1）长航班统计
flights %>%
  mutate(long_flight = (air_time >= 6 * 60)) %>%
  View()

flights %>%
  mutate(long_flight = (air_time >= 6 * 60)) %>%
  count(long_flight)

# 或者
flights %>%
  count(long_flight = air_time >= 6 * 60)

# 2）通过分组生成新的变量集
flights %>%
  group_by(date = make_date(year, month, day)) %>%
  summarise(flights_n = n(),
            air_time_mean = mean(air_time, na.rm = TRUE)) %>%
  ungroup()

# 3）随机获取样例
flights %>%
  slice_sample(n = 15)

flights %>%
  slice_sample(prop = 0.15)

# 4）创建日期列
flights %>%
  select(year, month, day) %>%
  mutate(date = make_date(year, month, day))

# 5）数字解析
numbers_1 <- tibble(number = c("#1", "Number8", "How are you 3"))
numbers_1 %>% mutate(number = parse_number(number))

# 6）多条件生成变量，按着循序执行和获取对应结果
flights %>%
  mutate(
    origin = case_when(
      (origin == "EWR") &
        dep_delay > 20 ~ "Newark International Airport - DELAYED",
      (origin == "EWR") &
        dep_delay <= 20 ~ "Newark International Airport - ON TIME DEPARTURE",
    )
  ) %>%
  count(origin)


# 7）一次性替换匹配到的所有模式
flights %>%
  mutate(origin = str_replace_all(
    origin,
    c("^EWR$" = "Newark International",    "^JFK$" = "John F. Kennedy International")
  )) %>%
  count(origin)

# 8）过滤数据集，行选择
flights_top_carriers <- flights %>%
  group_by(carrier) %>%
  filter(n() >= 10000) %>%
  ungroup()

flights_top_carriers %>% View

# 9) 抽取行，字符串检测
beginning_with_am <- airlines %>%
  filter(name %>% str_detect("^Am"))
beginning_with_am %>% View

# 10）补集和交集运算
data1 <- flights %>%
  anti_join(beginning_with_am, by = "carrier") %>% 
  View
data2 <- flights %>%
  inner_join(beginning_with_am, by = "carrier") %>% 
  View
nrow(data1)
nrow(data2)
nrow(flights)

# 11) fct_reorder
airline_names <- flights %>%
  left_join(airlines, by = "carrier")

airline_names %>%
  count(name) %>%
  ggplot(aes(name, n)) +
  geom_col()

airline_names %>%
  count(name) %>%
  mutate(name = fct_reorder(name, n)) %>%
  ggplot(aes(name, n)) +
  geom_col()

# 12) 坐标轴转换
airline_names %>%   
  count(name) %>%   
  mutate(name = fct_reorder(name, n)) %>%   
  ggplot(aes(name, n)) +   
  geom_col() +   
  coord_flip() 

# 13）Crossing操作
# 生成笛卡尔积
crossing(
  customer_channel = c("Bus", "Car"),
  customer_status = c("New", "Repeat"),
  spend_range = c("$0-$10", "$10-$20", "$20-$50", "$50+"))

# 14）基于自定义函数做分组
summary1 <- function(data, col_names, na.rm = TRUE) {
  data %>%
    summarise(across({{ col_names }},
                     list(
                       min = min,
                       max = max,
                       median = median,
                       mean = mean
                     ),
                     na.rm = na.rm,
                     .names = "{col}_{fn}"
    ))
}

airline_names %>%
  summary1(c(air_time, arr_delay)) %>% 
  View
airline_names %>%
  group_by(carrier) %>%
  summary1(c(air_time, arr_delay)) %>% 
  View

rm(list = ls())
# 数据集
set.seed(1045)
data <- data.frame(group = rep(c("P1", "P2", "P3"), each = 40),
                   values = c(rnorm(40, 0, 3), rnorm (40, 0, 6), rnorm (40, 1, 5)))
head(data)

# 方差分析模型
model <- aov(values ~ group, data = data)
summary(model)

# 做Tukey HSD 检验
TukeyHSD(model, conf.level = .95)

# 可视化分析
plot(TukeyHSD(model, conf.level=.95), las = 2)

data <- data.frame(x = c(1, 2, 5, 3, 5, 5, 9, 10, 12),
                   y = c(18, 10, 10, 20, 22, 13, 15, 16, 17))
data

# 1）添加线性趋势和置信区域
ggplot(data, aes(x = x, y = y)) +
  geom_point() +
  geom_smooth(method = lm, level = 0.99)

# 2）添加线性趋势和无置信区域
ggplot(data, aes(x = x, y = y)) +
  geom_point() +
  geom_smooth(
    method = lm,
    se = FALSE,
    col = 'blue',
    size = 2
  )

# 3）默认loess曲线
ggplot(data, aes(x = x, y = y)) +
  geom_point() +
  geom_smooth(se = FALSE)

graphics.off()
rm(list = ls())
options(warn = -1)
options(scipen = 999)
options(digits = 3)

# 创建数据框
manager <- c(1, 2, 3, 4, 5)
date <- c("10/24/08", "10/28/08", "10/1/08", "10/12/08", "5/1/09")
country <- c("US", "US", "UK", "UK", "UK")
gender <- c("M", "F", "F", "M", "F")
age <- c(32, 45, 25, 39, 99)
q1 <- c(5, 3, 3, 3, 2)
q2 <- c(4, 5, 5, 3, 2)
q3 <- c(5, 2, 5, 4, 1)
q4 <- c(5, 5, 5, NA, 2)
q5 <- c(5, 5, 2, NA, 1)
leadership <- data.frame(manager, date, country, gender, age,
                         q1, q2, q3, q4, q5, stringsAsFactors=FALSE)

# subset函数获取数据子集
# 1）选择所有age>=35或者age<24的行，保留变量q1, q3和q5
new_df1 <- subset(leadership, 
                 age >= 35 | age < 24, 
                 select = c(q1, q2, q5))
new_df1

# 2）选择所有25岁以上的男性，保留gender到q5的变量集
new_df2 <- subset(leadership,
                  gender == 'M' & age > 25,
                  select = gender:q5)
new_df2

# 第2集 研究the office数据集，情感分析和词云图
# 工作空间管理和配置
graphics.off()
rm(list = ls)
options(warn = -1)
options(scipen = 999)
options(digits = 3)

# R包
library(tidyr)
library(stringr)
library(schrute)
library(tidytext)
library(wordcloud)
library(ggplot2)
library(dplyr)
library(reshape2)

# 数据导入
office_ratings <-
  readr::read_csv('./data/tidytuesday/data/2020/2020-03-17/office_ratings.csv')

schrute <- schrute::theoffice
head(schrute) %>% View

# 查看数据
dplyr::glimpse(schrute)
dplyr::glimpse(office_ratings)

# 准备数据
token.schrute <-  schrute %>%
  tidytext::unnest_tokens(word, text)
dplyr::glimpse(token.schrute)

# Remove stop words
stop_words <-  tidytext::stop_words
tidy.token.schrute <-  token.schrute %>%
  dplyr::anti_join(stop_words, by = 'word')

# Most common words
tidy.token.schrute %>% # 169,835 observations
  dplyr::count(word, sort = TRUE)


# 词频分析
p1 <-  tidy.token.schrute %>%
  dplyr::count(word, sort = TRUE) %>%
  dplyr::filter(n > 400) %>%
  dplyr::mutate(word = stats::reorder(word, n)) %>%
  ggplot2::ggplot(ggplot2::aes(word, n)) +
  ggplot2::geom_col() +
  ggplot2::xlab(NULL) +
  ggplot2::coord_flip() +
  ggplot2::theme_minimal()

ggsave(filename = './figs/Most_used_words.pdf',
       plot = p1,
       scale = .6)

# 情感分析
# 导入情感词典
sentiments <-
  get_sentiments("bing") # Codes words as positive or negative (Bing Liu). NA for neutral.
dplyr::glimpse(sentiments)
unique(sentiments$sentiment)

# Sentiment by season
schrute.sentiment <-  tidy.token.schrute %>%
  dplyr::left_join(sentiments) %>%
  dplyr::count(episode_name, sentiment) %>%
  spread(sentiment, n, fill = 0) %>% # fill missing values w/ 0
  mutate(sentimentc = positive - negative) %>% # pos value means more words had positive connatation than neg
  dplyr::select(episode_name, sentimentc, negative, positive, neutral =
                  `<NA>`)

bing_word_counts <-  tidy.token.schrute %>%
  inner_join(sentiments %>% filter(sentiment == 'positive' |
                                     sentiment == 'negative')) %>%
  count(word, sentiment, sort = TRUE)
bing_word_counts

p2 <- bing_word_counts %>%
  filter(n > 150) %>%
  mutate(n = ifelse(sentiment == 'negative', -n, n)) %>%
  mutate(word = reorder(word, n)) %>%
  ggplot(aes(word, n, fill = sentiment)) +
  geom_col() +
  coord_flip() +
  labs(y = 'Contribution to sentiment analysis', x = 'Word') +
  theme_bw() +
  theme(
    legend.position = 'none',
    axis.text.x = element_text(size = 12),
    axis.title.x = element_text(size = 14, face = 'bold'),
    axis.text.y = element_text(size = 12),
    axis.title.y = element_text(size = 14, face = 'bold')
  )
p2

# word cloud
pdf('./figs/comparison_cloud.pdf',
    width = 4,
    height = 4)
bing_word_counts %>%
  acast(word ~ sentiment, value.var = 'n', fill = 0) %>%
  comparison.cloud(colors = c("#F8766D", "#00BFC4"),
                   max.words = 100)
dev.off()


# ?Is there a relationship between word sentiment and episode rating?
office_ratings <-  office_ratings %>%
  dplyr::select(season, episode, episode_name = title, imdb_rating)

sent.rating <-  schrute.sentiment %>%
  inner_join(office_ratings, by = 'episode_name') %>%
  mutate(season = as.factor(season))

ggplot(data = sent.rating, aes(x = sentimentc, y = imdb_rating, color = season)) +
  geom_point()

sent.rating %>%
  group_by(season) %>%
  summarize(rating_ave = mean(imdb_rating),
            sentiment_ave = mean(sentimentc)) %>%
  ggplot(data = ., aes(x = sentiment_ave, y = rating_ave, color = season)) +
  geom_point()

# Descriptive fig of pos & neg characters
glimpse(tidy.token.schrute)

char.sentiment <-  tidy.token.schrute %>%
  inner_join(sentiments, by = 'word') %>%
  count(character, sentiment) %>%
  spread(sentiment, n, fill = 0) %>%
  mutate(sentimentc = positive - negative)

p3 <-  char.sentiment %>%
  filter(negative + positive > 300) %>%
  mutate(sent_dummy = ifelse(sentimentc < 0, 'More Negative', 'More Positive')) %>%
  mutate(character = reorder(character, sentimentc)) %>%
  ggplot(aes(character, sentimentc, fill = sent_dummy)) +
  geom_col() +
  coord_flip() +
  labs(y = 'Emotional Charge of Dialogue n (Positive - Negative Words)', x = 'Character') +
  theme_bw() +
  theme(
    legend.position = 'none',
    axis.text.x = element_text(size = 12),
    axis.title.x = element_text(size = 14, face = 'bold'),
    axis.text.y = element_text(size = 12),
    axis.title.y = element_text(size = 14, face = 'bold')
  )

ggsave(
  filename = './figs/word_bar.pdf',
  plot = p2,
  width = 4.5,
  height = 4
)
ggsave(
  filename = './figs/characters_sentiment.pdf',
  plot = p3,
  width = 9,
  height = 6
)