Цель: Исследовать данные о поведении пользователей онлайн-кинотеатра KION с использованием статистических методов.
Данные представляют собой журнал взаимодействий пользователей с контентом за период с 13 марта 2021 года по 22 августа 2022 года. Для анализа использованы следующие переменные:
age – возрастная группа пользователя (категориальная, 6 групп).
income – уровень дохода пользователя (категориальная, 6 групп).
total_dur – общая продолжительность всех просмотров в секундах (количественная).
rating_kp – рейтинг фильма на Кинопоиске (количественная, от 0 до 10).
sex – пол пользователя (категориальная, М/Ж).
library(dplyr)
##
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
library(readxl)
library(ggplot2)
library(car)
## Loading required package: carData
##
## Attaching package: 'car'
## The following object is masked from 'package:dplyr':
##
## recode
library(moments)
library(knitr)
library(vcd)
## Loading required package: grid
library(viridis)
## Loading required package: viridisLite
library(scales)
##
## Attaching package: 'scales'
## The following object is masked from 'package:viridis':
##
## viridis_pal
df <- read_excel('/Users/privet/Downloads/task_data.xlsx')
records <- nrow(df)
cat("Количество записей:", records, "\n")
## Количество записей: 15768
unique_users <- length(unique(df$user_id))
cat("Количество уникальных пользователей:", unique_users, "\n")
## Количество уникальных пользователей: 15238
unique_items <- length(unique(df$item_id))
cat("Количество уникальных единиц контента (по item_id):", unique_items, "\n")
## Количество уникальных единиц контента (по item_id): 3204
films <- length(unique(df$item_id[df$content_type == "film"]))
series <- length(unique(df$item_id[df$content_type == "series"]))
cat("Количество уникальных фильмов (по item_id):", films, "\n")
## Количество уникальных фильмов (по item_id): 2644
cat("Количество уникальных сериалов (по item_id):", series, "\n")
## Количество уникальных сериалов (по item_id): 560
3.1 Возраст (age)
Возраст представлен в виде шести групп.
Распределение по группам:
age_18_24 – 18-24 года
age_25_34 – 25-34 года
age_35_44 – 35-44 года
age_45_54 – 45-54 года
age_55_64 – 55-64 года
age_65_inf – 65 лет и старше
df %>%
filter(!is.na(age)) %>%
count(age) %>%
mutate(age = case_when(
age == "age_18_24" ~ "18-24 лет",
age == "age_25_34" ~ "25-34 лет",
age == "age_35_44" ~ "35-44 лет",
age == "age_45_54" ~ "45-54 лет",
age == "age_55_64" ~ "55-64 лет",
age == "age_65_inf" ~ "65+ лет",
TRUE ~ age
)) %>%
ggplot(aes(x=age, y=n)) +
geom_col(fill="#3498db", alpha=0.8, width=0.7) +
geom_text(aes(label=n), vjust=-0.5, size=3.5, color="gray30") +
theme_minimal(base_size=12) +
labs(title="Количество просмотров по возрастным группам",
subtitle="Распределение просмотров среди возрастных категорий",
x="Возрастная группа",
y="Число просмотров") +
theme(
plot.title = element_text(face="bold", size=16, hjust=0.5),
plot.subtitle = element_text(color="gray40", hjust=0.5, size=10),
axis.title = element_text(face="bold"),
axis.text = element_text(color="gray30"),
panel.grid.major.x = element_blank(),
panel.grid.minor = element_blank()
) +
scale_y_continuous(expand = expansion(mult = c(0, 0.1)))
Вывод: Наибольшее количество просмотров приходится на возрастные группы 25-34 года и 35-44 года. Наименьшее – на группу 65+. Это может свидетельствовать о том, что основная аудитория KION – люди среднего возраста.
3.2 Доход (income)
Доход также разбит на группы (в рублях):
income_0_20 – 0 – 20,000
income_20_40 – 20,000 – 40,000
income_40_60 – 40,000 – 60,000
income_60_90 – 60,000 – 90,000
income_90_150 – 90,000 – 150,000
income_150_inf – более 150,000
График: столбчатая диаграмма частот доходных групп.
df %>%
filter(!is.na(income)) %>%
count(income) %>%
mutate(
income = case_when(
income == "income_0_20" ~ "0-20k",
income == "income_20_40" ~ "20-40k",
income == "income_40_60" ~ "40-60k",
income == "income_60_90" ~ "60-90k",
income == "income_90_150" ~ "90-150k",
income == "income_150_inf" ~ "150k+",
TRUE ~ income
),
# Создаем фактор с правильным порядком
income = factor(income, levels = c("0-20k", "20-40k", "40-60k",
"60-90k", "90-150k", "150k+"))
) %>%
ggplot(aes(x=income, y=n)) +
geom_col(fill="#e74c3c", alpha=0.8, width=0.7) +
geom_text(aes(label=n), vjust=-0.5, size=3.5, color="gray30") +
theme_minimal(base_size=12) +
labs(title="Распределение просмотров по уровню дохода",
subtitle="Количество просмотров в различных доходных группах",
x="Доходная группа (тыс. $)",
y="Число просмотров") +
theme(
plot.title = element_text(face="bold", size=16, hjust=0.5),
plot.subtitle = element_text(color="gray40", hjust=0.5, size=10),
axis.title = element_text(face="bold"),
axis.text = element_text(color="gray30"),
axis.text.x = element_text(angle=45, hjust=1),
panel.grid.major.x = element_blank(),
panel.grid.minor = element_blank()
) +
scale_y_continuous(expand = expansion(mult = c(0, 0.1)))
Вывод: Преобладают пользователи с доходом 20 000 – 40 000 руб. и 40 000 – 60 000 руб., что соответствует среднему классу. Меньше всего просмотров у группы с доходом более 150 000 руб.
3.3 Пол (sex)
sex_clean <- df$sex[!is.na(df$sex) & df$sex != ""]
sex_table <- table(sex_clean)
# Процентное соотношение
sex_percent <- prop.table(sex_table) * 100
# Красивый вывод
cat("\nПроцентное соотношение мужчин и женщин:\n")
##
## Процентное соотношение мужчин и женщин:
cat("Мужчины:", round(sex_percent["М"], 2), "%\n")
## Мужчины: 51.58 %
cat("Женщины:", round(sex_percent["Ж"], 2), "%\n")
## Женщины: 48.42 %
Вывод: Доля мужчин незначительно выше.
3.4 Продолжительность просмотров (total_dur)
# Извлекаем продолжительность просмотров
total_dur <- df$total_dur
# Основные статистики
min_dur <- min(total_dur, na.rm = TRUE)
max_dur <- max(total_dur, na.rm = TRUE)
mean_dur <- mean(total_dur, na.rm = TRUE)
median_dur <- median(total_dur, na.rm = TRUE)
# Перевод в удобный формат (часы:минуты:секунды)
seconds_to_time <- function(seconds) {
hours <- floor(seconds / 3600)
minutes <- floor((seconds %% 3600) / 60)
secs <- round(seconds %% 60)
return(paste(hours, "ч", minutes, "мин", secs, "сек"))
}
# Вывод результатов
cat("Минимальная длительность:", min_dur, "сек (", seconds_to_time(min_dur), ")\n")
## Минимальная длительность: 1 сек ( 0 ч 0 мин 1 сек )
cat("Максимальная длительность:", max_dur, "сек (", seconds_to_time(max_dur), ")\n")
## Максимальная длительность: 3502510 сек ( 972 ч 55 мин 10 сек )
cat("Среднее время просмотра(арифметическое):", round(mean_dur), "сек (", seconds_to_time(mean_dur), ") - сильно зависит от количества очень длительных просмотров, которые завышают среднее значениe)\n")
## Среднее время просмотра(арифметическое): 8558 сек ( 2 ч 22 мин 38 сек ) - сильно зависит от количества очень длительных просмотров, которые завышают среднее значениe)
cat("Медиана:", median_dur, "сек (", seconds_to_time(median_dur),")", "- устойчива к выбросам и показывает, что половина всех просмотров длится менее 48 минут\n\n")
## Медиана: 2882 сек ( 0 ч 48 мин 2 сек ) - устойчива к выбросам и показывает, что половина всех просмотров длится менее 48 минут
# Дополнительная информация
cat("Дополнительные статистики\n")
## Дополнительные статистики
cat("Просмотров < 10 сек:", sum(total_dur < 10, na.rm = TRUE), "\n")
## Просмотров < 10 сек: 525
cat("Просмотров > 2 часов (>7200 сек):", sum(total_dur > 7200, na.rm = TRUE), "\n")
## Просмотров > 2 часов (>7200 сек): 3887
cat("Просмотров > 1 дня (>86400 сек):", sum(total_dur > 86400, na.rm = TRUE), "\n")
## Просмотров > 1 дня (>86400 сек): 174
# Квартили
quantiles <- quantile(total_dur, probs = c(0.25, 0.5, 0.75), na.rm = TRUE)
cat("\nКвартили:\n")
##
## Квартили:
cat("25% просмотров короче:", quantiles[1], "сек (", seconds_to_time(quantiles[1]), ")\n")
## 25% просмотров короче: 379.75 сек ( 0 ч 6 мин 20 сек )
cat("50% просмотров короче:", quantiles[2], "сек (", seconds_to_time(quantiles[2]), ")\n")
## 50% просмотров короче: 2882 сек ( 0 ч 48 мин 2 сек )
cat("75% просмотров короче:", quantiles[3], "сек (", seconds_to_time(quantiles[3]), ")\n")
## 75% просмотров короче: 7134.25 сек ( 1 ч 58 мин 54 сек )
Я использую логарифмическую шкалу на графике, поскольку переменная total_dur имеет сильную правостороннюю асимметрию: огромное количество коротких просмотров (несколько минут) соседствует с небольшим числом очень длинных (часы и дни). В таких случаях линейная шкала делает график нечитаемым — длинные просмотры “растягивают” ось, а основная масса данных сжимается в левой части. Логарифмирование преобразует шкалу, делая распределение более симметричным и позволяя визуально сравнить как короткие, так и длинные просмотры.
# Сначала создадим датафрейм с нужными данными
df_duration <- df %>%
filter(total_dur > 0) %>%
mutate(log_duration = log10(total_dur))
# Вычислим среднее и медиану в логарифмической шкале
mean_log <- mean(log10(df$total_dur[df$total_dur > 0]), na.rm = TRUE)
median_log <- median(log10(df$total_dur[df$total_dur > 0]), na.rm = TRUE)
# Построим гистограмму
df_duration %>%
ggplot(aes(x = log_duration)) +
geom_histogram(fill = "#87CEEB", alpha = 0.8, color = "white", bins = 30) +
theme_minimal(base_size = 12) +
labs(title = "Распределение продолжительности просмотров",
subtitle = "Гистограмма длительности в логарифмической шкале",
x = "log10(длительность в секундах)",
y = "Частота") +
theme(
plot.title = element_text(face = "bold", size = 16, hjust = 0.5),
plot.subtitle = element_text(color = "gray40", hjust = 0.5, size = 10),
axis.title = element_text(face = "bold"),
axis.text = element_text(color = "gray30"),
panel.grid.minor = element_blank()
) +
scale_x_continuous(breaks = seq(0, 6, by = 1)) +
scale_y_continuous(limits = c(0, 4000), expand = c(0, 0))
Выводы из распределения:
На гистограмме четко выделяются три основных паттерна поведения:
Доминирование коротких просмотров (log10 < 3, т.е. < 1000 секунд). Это самая многочисленная группа. Пользователи часто ограничиваются знакомством с началом фильма или просмотром одной-двух серий сериала. Это может указывать на “пробное” потребление или на то, что контент не смог удержать внимание зрителя.
Просмотры средней длительности (log10 3.5–4, т.е. ~30–100 минут). Эта группа соответствует просмотру короткометражных фильмов или значительной части полнометражного фильма. Здесь представлены зрители, которые уже достаточно вовлечены, но по разным причинам не досматривают до конца.
Длинные просмотры (log10 > 4.5, т.е. > 3-4 часов). Наименее многочисленная категория. Эти значения, скорее всего, относятся к полнометражным фильмам, просмотренным целиком, или к нескольким сериям одного сериала подряд (“binge-watching”).
3.5 Рейтинг Кинопоиска (rating_kp)
rating <- df$rating_kp[!is.na(df$rating_kp)]
cat("Минимальный рейтинг:", min(rating), "\n")
## Минимальный рейтинг: 0
cat("Максимальный рейтинг:", max(rating), "\n")
## Максимальный рейтинг: 9.2
cat("Средний рейтинг:", round(mean(rating), 2), "\n")
## Средний рейтинг: 6.62
cat("Медианный рейтинг:", median(rating), "\n")
## Медианный рейтинг: 6.8
cat("Стандартное отклонение:", round(sd(rating), 2), "\n")
## Стандартное отклонение: 1.23
# Квартили
q <- quantile(rating, probs = c(0.25, 0.5, 0.75))
cat("\nКвартили:\n")
##
## Квартили:
cat("25% фильмов имеют рейтинг ниже:", q[1], "\n")
## 25% фильмов имеют рейтинг ниже: 6.1
cat("50% фильмов имеют рейтинг ниже (медиана):", q[2], "\n")
## 50% фильмов имеют рейтинг ниже (медиана): 6.8
cat("75% фильмов имеют рейтинг ниже:", q[3], "\n")
## 75% фильмов имеют рейтинг ниже: 7.5
rating_groups <- cut(rating,
breaks = c(0, 4, 7, 10),
labels = c("0-4", "5-7", "8-10"),
include.lowest = TRUE,
right = TRUE)
# сколько фильмов в каждой группе
group_counts <- table(rating_groups)
group_percents <- prop.table(group_counts) * 100
# таблица результатов
result_table <- data.frame(
Рейтинг = names(group_counts),
Количество = as.vector(group_counts),
Процент = round(as.vector(group_percents), 1)
)
print("Распределение фильмов по рейтингу:")
## [1] "Распределение фильмов по рейтингу:"
print(result_table)
## Рейтинг Количество Процент
## 1 0-4 276 2.2
## 2 5-7 7725 61.1
## 3 8-10 4636 36.7
# Рассчитываем статистики на всех данных (кроме NA)
mean_rating <- mean(df$rating_kp, na.rm = TRUE)
median_rating <- median(df$rating_kp, na.rm = TRUE)
# Строим график с coord_cartesian вместо limits
df %>%
filter(!is.na(rating_kp)) %>%
ggplot(aes(x = rating_kp)) +
geom_histogram(
fill = "#2ecc71",
alpha = 0.8,
color = "white",
bins = 15
) +
geom_vline(
xintercept = mean_rating,
color = "#e74c3c",
linetype = "dashed",
size = 1
) +
geom_vline(
xintercept = median_rating,
color = "#f39c12",
linetype = "dashed",
size = 1
) +
annotate(
"text",
x = mean_rating,
y = Inf,
label = paste("Среднее:", round(mean_rating, 2)),
vjust = 2,
hjust = -0.1,
color = "black",
size = 3.5
) +
annotate(
"text",
x = median_rating,
y = Inf,
label = paste("Медиана:", round(median_rating, 2)),
vjust = 4,
hjust = -0.1,
color = "black",
size = 3.5
) +
theme_minimal(base_size = 12) +
labs(
title = "Распределение рейтинга Кинопоиска",
subtitle = "Гистограмма рейтингов фильмов",
x = "Рейтинг",
y = "Количество фильмов"
) +
theme(
plot.title = element_text(face = "bold", size = 16, hjust = 0.5),
plot.subtitle = element_text(color = "gray40", hjust = 0.5, size = 10),
axis.title = element_text(face = "bold"),
axis.text = element_text(color = "gray30"),
panel.grid.minor = element_blank()
) +
coord_cartesian(xlim = c(0, 10)) + # обрезает видимую область, не удаляя данные
scale_x_continuous(breaks = seq(0, 10, by = 2)) +
scale_y_continuous(expand = expansion(mult = c(0, 0.1)))
## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ℹ Please use `linewidth` instead.
## This warning is displayed once per session.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
Вывод: Пользователи отдают преимущество контенту со средним и выше среднего рейтингом (6-7). Фильмы с очень низким рейтингом (<4) встречаются крайне редко.
Нулевая гипотеза (H₀): Нет статистически значимых различий в продолжительности просмотра контента между группами пользователей с разным уровнем дохода.
Альтернативная гипотеза (H₁): Существуют статистически значимые различия в продолжительности просмотра контента между группами пользователей с разным уровнем дохода.
# Подготовка данных
df_clean <- df %>%
filter(!is.na(income) & income != "" & !is.na(total_dur) & total_dur > 0)
# income в фактор
df_clean$income <- as.factor(df_clean$income)
cat("ПЕРВИЧНЫЙ ОСМОТР ДАННЫХ\n")
## ПЕРВИЧНЫЙ ОСМОТР ДАННЫХ
cat("Всего наблюдений:", nrow(df_clean), "\n")
## Всего наблюдений: 12515
cat("Диапазон длительности: от", min(df_clean$total_dur)/60, "мин до",
max(df_clean$total_dur)/3600, "часов\n")
## Диапазон длительности: от 0.01666667 мин до 972.9194 часов
cat("Среднее:", mean(df_clean$total_dur)/60, "минут\n")
## Среднее: 134.304 минут
cat("Медиана:", median(df_clean$total_dur)/60, "минут\n\n")
## Медиана: 45.41667 минут
поскольку в данных есть выбросы, нужно от них избавиться.
# Удаление выбросов по методу межквартильного размаха (IQR)
remove_outliers <- function(df, group_col, value_col) {
df %>%
group_by(!!sym(group_col)) %>%
mutate(
Q1 = quantile(!!sym(value_col), 0.25, na.rm = TRUE),
Q3 = quantile(!!sym(value_col), 0.75, na.rm = TRUE),
IQR = Q3 - Q1,
lower_bound = Q1 - 1.5 * IQR,
upper_bound = Q3 + 1.5 * IQR
) %>%
filter(!!sym(value_col) >= lower_bound & !!sym(value_col) <= upper_bound) %>%
select(-Q1, -Q3, -IQR, -lower_bound, -upper_bound)
}
df_no_outliers <- remove_outliers(df_clean, "income", "total_dur")
cat("Было наблюдений:", nrow(df_clean), "\n")
## Было наблюдений: 12515
cat("Стало наблюдений:", nrow(df_no_outliers), "\n")
## Стало наблюдений: 11489
cat("Удалено выбросов:", nrow(df_clean) - nrow(df_no_outliers),
"(", round((nrow(df_clean) - nrow(df_no_outliers))/nrow(df_clean)*100, 1), "%)\n\n")
## Удалено выбросов: 1026 ( 8.2 %)
cat("После удаления выбросов:\n")
## После удаления выбросов:
cat("Среднее:", mean(df_no_outliers$total_dur)/60, "минут\n")
## Среднее: 58.89118 минут
cat("Медиана:", median(df_no_outliers$total_dur)/60, "минут\n")
## Медиана: 33.85 минут
cat("Максимум:", max(df_no_outliers$total_dur)/60, "минут\n\n")
## Максимум: 342.0667 минут
# 4.1 Визуальная проверка нормальности (Q-Q plots)
par(mfrow = c(2, 3))
for(g in levels(df_no_outliers$income)) {
data_subset <- df_no_outliers$total_dur[df_no_outliers$income == g]
if(length(data_subset) > 0) {
qqnorm(data_subset, main = paste("Q-Q plot для", g))
qqline(data_subset, col = "red")
}
}
par(mfrow = c(1, 1))
# 4.2 Тест Ливиня
levene_result <- leveneTest(total_dur ~ income, data = df_no_outliers)
print(levene_result)
## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 5 5.6556 3.257e-05 ***
## 11483
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Дисперсии неоднородны (p < 0.05), поэтому ANOVA может быть ненадежна. Используем тест Крускала-Уоллиса.
kruskal_result <- kruskal.test(total_dur ~ income, data = df_no_outliers)
print(kruskal_result)
##
## Kruskal-Wallis rank sum test
##
## data: total_dur by income
## Kruskal-Wallis chi-squared = 2.8582, df = 5, p-value = 0.7218
p_value <- kruskal_result$p.value
test_used <- "Крускала-Уоллиса"
ИНТЕРПРЕТАЦИЯ РЕЗУЛЬТАТОВ ТЕСТА
Тест Крускала-Уоллиса:
Хи-квадрат = 2.8582
p-значение = 0.7218 > 0.05
- НЕТ статистически значимых различий между группами
- Нулевая гипотеза НЕ ОТВЕРГАЕТСЯ
- Доход НЕ влияет на продолжительность просмотра контента
df_summary <- df_no_outliers %>%
group_by(income) %>%
summarise(
count = n(),
mean_duration = mean(total_dur),
mean_minutes = round(mean(total_dur)/60, 1),
median_duration = median(total_dur),
median_minutes = round(median(total_dur)/60, 1),
sd_duration = sd(total_dur),
min_minutes = round(min(total_dur)/60, 1),
max_minutes = round(max(total_dur)/60, 1)
) %>%
arrange(desc(median_duration))
print(df_summary)
## # A tibble: 6 × 9
## income count mean_duration mean_minutes median_duration median_minutes
## <fct> <int> <dbl> <dbl> <dbl> <dbl>
## 1 income_40_60 3664 3558. 59.3 2096. 34.9
## 2 income_90_150 217 4183. 69.7 2048 34.1
## 3 income_20_40 6357 3484. 58.1 2031 33.9
## 4 income_60_90 980 3619. 60.3 1892 31.5
## 5 income_0_20 258 3411. 56.9 1669 27.8
## 6 income_150_inf 13 6071. 101. 613 10.2
## # ℹ 3 more variables: sd_duration <dbl>, min_minutes <dbl>, max_minutes <dbl>
median_values <- df_summary$median_minutes
names(median_values) <- df_summary$income
# группы с наибольшими и наименьшими значениями
max_group <- df_summary[which.max(df_summary$median_duration), ]
min_group <- df_summary[which.min(df_summary$median_duration), ]
cat("Наибольшая медианная продолжительность:\n")
## Наибольшая медианная продолжительность:
cat(" Группа:", as.character(max_group$income), "\n")
## Группа: income_40_60
cat(" Медиана:", max_group$median_minutes, "минут\n")
## Медиана: 34.9 минут
cat(" Среднее:", max_group$mean_minutes, "минут\n")
## Среднее: 59.3 минут
cat(" Наблюдений:", max_group$count, "\n\n")
## Наблюдений: 3664
cat("Наименьшая медианная продолжительность:\n")
## Наименьшая медианная продолжительность:
cat(" Группа:", as.character(min_group$income), "\n")
## Группа: income_150_inf
cat(" Медиана:", min_group$median_minutes, "минут\n")
## Медиана: 10.2 минут
cat(" Среднее:", min_group$mean_minutes, "минут\n")
## Среднее: 101.2 минут
cat(" Наблюдений:", min_group$count, "\n\n")
## Наблюдений: 13
# Разница между группами
diff_minutes <- max_group$median_minutes - min_group$median_minutes
diff_percent <- round((max_group$median_duration - min_group$median_duration) /
min_group$median_duration * 100, 1)
cat("Разница между группами:", diff_minutes, "минут", "\n")
## Разница между группами: 24.7 минут
Но эта разница статистически НЕ значима (p = 0.7218)
# Сначала создадим фактор с правильным порядком уровней
desired_order <- c("0-20" = "income_0_20",
"20-40" = "income_20_40",
"40-60" = "income_40_60",
"60-90" = "income_60_90",
"90-150" = "income_90_150",
"150+" = "income_150_inf")
# Создаем новый столбец с красивыми названиями и правильным порядком
df_plot <- df_no_outliers %>%
mutate(
income_pretty = case_when(
income == "income_0_20" ~ "0-20 тыс.",
income == "income_20_40" ~ "20-40 тыс.",
income == "income_40_60" ~ "40-60 тыс.",
income == "income_60_90" ~ "60-90 тыс.",
income == "income_90_150" ~ "90-150 тыс.",
income == "income_150_inf" ~ "150+ тыс.",
TRUE ~ income
),
income_pretty = factor(income_pretty,
levels = c("0-20 тыс.", "20-40 тыс.", "40-60 тыс.",
"60-90 тыс.", "90-150 тыс.", "150+ тыс."))
)
# Создаем улучшенный график с фиксированной шириной столбцов
p1 <- ggplot(df_plot, aes(x = income_pretty, y = total_dur/60, fill = income_pretty)) +
geom_boxplot(
alpha = 0.8,
width = 0.7, # Фиксированная ширина столбцов
outlier.color = "darkred",
outlier.shape = 16,
outlier.size = 2,
outlier.alpha = 0.5,
notch = FALSE,
varwidth = FALSE # Отключаем переменную ширину
) +
stat_summary(
fun = mean,
geom = "point",
shape = 18,
size = 4,
color = "darkblue",
fill = "darkblue"
) +
scale_fill_brewer(palette = "Blues") +
theme_minimal(base_size = 12) +
labs(
title = "Продолжительность просмотра по группам дохода",
subtitle = paste("Тест Крускала-Уоллиса: p = 0.7218 (статистически не значимо)"),
x = "Доходная группа (тыс. руб.)",
y = "Продолжительность просмотра (минуты)",
caption = "Красная пунктирная линия - общая медиана\nСиние ромбы - средние значения"
) +
theme(
axis.text.x = element_text(angle = 0, hjust = 0.5, face = "bold"),
axis.title = element_text(face = "bold"),
plot.title = element_text(face = "bold", size = 14),
plot.subtitle = element_text(color = "darkred", size = 11),
plot.caption = element_text(color = "gray50", size = 9, hjust = 0, face = "plain"),
legend.position = "none",
panel.grid.minor = element_blank(),
panel.grid.major.x = element_blank(),
panel.border = element_blank(),
axis.line = element_line(color = "gray70"),
plot.background = element_rect(fill = "white", color = NA)
) +
# Добавляем только линию общей медианы, без подписи
geom_hline(
yintercept = median(df_no_outliers$total_dur)/60,
linetype = "dashed",
color = "red",
alpha = 0.9,
size = 1
)
print(p1)
# Сначала подготовим данные с правильным порядком групп
df_summary <- df_no_outliers %>%
group_by(income) %>%
summarise(median_duration = median(total_dur, na.rm = TRUE)) %>%
mutate(
median_minutes = round(median_duration/60, 1),
income_pretty = case_when(
income == "income_0_20" ~ "0-20 тыс.",
income == "income_20_40" ~ "20-40 тыс.",
income == "income_40_60" ~ "40-60 тыс.",
income == "income_60_90" ~ "60-90 тыс.",
income == "income_90_150" ~ "90-150 тыс.",
income == "income_150_inf" ~ "150+ тыс.",
TRUE ~ income
),
income_pretty = factor(income_pretty,
levels = c("0-20 тыс.", "20-40 тыс.", "40-60 тыс.",
"60-90 тыс.", "90-150 тыс.", "150+ тыс."))
) %>%
arrange(income_pretty) # Сортируем по порядку
# Создаем столбчатую диаграмму
p2 <- ggplot(df_summary, aes(x = income_pretty, y = median_minutes, fill = income_pretty)) +
geom_col(
width = 0.7, # Такая же ширина, как у боксплотов
alpha = 0.8, # Такая же прозрачность
color = "gray30", # Легкая обводка для аккуратности
linewidth = 0.3
) +
geom_text(
aes(label = paste0(median_minutes, " мин")),
vjust = -0.8, # Немного выше, чтобы не сливалось со столбцом
size = 4,
fontface = "bold",
color = "darkblue"
) +
scale_fill_brewer(palette = "Blues") + # Такая же цветовая гамма
theme_minimal(base_size = 12) +
labs(
title = "Медианная продол-ть просмотра по группам дохода",
subtitle = paste("Различия статистически НЕ значимы (p = 0.7218)"),
x = "Доходная группа (тыс. руб.)",
y = "Медианная продолжительность (минуты)",
caption = "Столбцы показывают медианные значения по группам"
) +
theme(
axis.text.x = element_text(angle = 0, hjust = 0.5, face = "bold"), # Горизонтальный текст
axis.title = element_text(face = "bold"),
plot.title = element_text(face = "bold", size = 14),
plot.subtitle = element_text(color = "darkred", size = 11),
plot.caption = element_text(color = "gray50", size = 9, hjust = 0, face = "plain"),
legend.position = "none",
panel.grid.minor = element_blank(),
panel.grid.major.x = element_blank(), # Убираем вертикальные линии сетки
panel.border = element_blank(),
axis.line = element_line(color = "gray70"),
plot.background = element_rect(fill = "white", color = NA)
) +
# Добавляем общую медиану для контекста
geom_hline(
yintercept = median(df_no_outliers$total_dur)/60,
linetype = "dashed",
color = "red",
alpha = 0.7,
size = 0.8
) +
scale_y_continuous(
limits = c(0, max(df_summary$median_minutes) * 1.15), # Немного места для подписей
expand = expansion(mult = c(0, 0.05)) # Убираем отступ снизу, оставляем сверху
)
print(p2)
НЕТ, статистически значимых различий не обнаружено (p = 0.72). Уровень дохода не влияет на то, как долго пользователи смотрят контент.
Формально можно выделить группу с наибольшим медианным временем просмотра, однако эти различия статистически не значимы, поэтому данный вывод ненадежен.
Аналогично, формально можно определить группу с наименьшим медианным временем, но эти различия статистически не значимы.
ОБЩИЙ ВЫВОД:
После удаления выбросов и применения робастного теста Крускала-Уоллиса не обнаружено статистически значимых различий в продолжительности просмотра между группами с разным уровнем дохода. Доход пользователя не является фактором, определяющим длительность просмотра контента в KION, что подтверждает нулевую гипотезу.
Нулевая гипотеза (H₀): Нет различий в средней продолжительности просмотра контента между мужчинами и женщинами.
Альтернативная гипотеза (H₁): Существуют статистически значимые различия в средней продолжительности просмотра контента между мужчинами и женщинами.
Для сравнения двух независимых групп (мужчины и женщины) я буду использовать двухвыборочный t-тест (t-test) или его непараметрический аналог U-тест Манна-Уитни (Wilcoxon rank-sum test).
У нас две независимые группы (мужчины и женщины) Зависимая переменная — количественная (продолжительность просмотра в секундах) t-тест — параметрический, требует нормальности распределения U-тест Манна-Уитни — непараметрический, не требует нормальности
Исходя из прошлого задания, распределение продолжительности просмотра скошено и содержит выбросы. Поэтому я буду использовать U-тест Манна-Уитни как более робастный, но для полноты проведу оба теста.
# убираем пропуски
df_clean <- df %>%
filter(!is.na(sex) & sex != "" & !is.na(total_dur) & total_dur > 0)
# Преобразуем пол в фактор
df_clean$sex <- as.factor(df_clean$sex)
cat("РАЗМЕР ГРУПП\n")
## РАЗМЕР ГРУПП
table(df_clean$sex)
##
## Ж М
## 6051 6446
cat("\n")
# Логарифмированная шкала для визуализации
boxplot(log10(total_dur) ~ sex, data = df_clean,
main = "Продолжительность просмотра (log10)",
ylab = "log10(секунды)", xlab = "Пол",
col = c("lightblue", "pink"))
par(mfrow = c(1, 1))
# Проверка нормальности (Q-Q plots)
par(mfrow = c(1, 2))
for(g in levels(df_clean$sex)) {
data_subset <- df_clean$total_dur[df_clean$sex == g]
qqnorm(data_subset, main = paste("Q-Q plot для пола", g))
qqline(data_subset, col = "red")
}
par(mfrow = c(1, 1))
# Проверка гомогенности дисперсий (тест Ливиня)
levene_result <- leveneTest(total_dur ~ sex, data = df_clean)
print(levene_result)
## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 1 1.8006 0.1797
## 12495
Дисперсии однородны (p = 0.1797 > 0.05), значит мы можем использовать t-тест (параметрический) как основной, так как одно из важных допущений выполнено. Но можно заметить, что в данных большое количество выбросов, от которых нужно избавиться.
# Функция для удаления ВЕРХНИХ выбросов
remove_upper_outliers <- function(df, group_col, value_col) {
df %>%
group_by(!!sym(group_col)) %>%
mutate(
Q3 = quantile(!!sym(value_col), 0.75, na.rm = TRUE),
IQR = IQR(!!sym(value_col), na.rm = TRUE),
upper_bound = Q3 + 1.5 * IQR # только верхняя граница
) %>%
filter(!!sym(value_col) <= upper_bound) %>% # оставляем все, что НИЖЕ верхней границы
select(-Q3, -IQR, -upper_bound)
}
df_no_outliers <- remove_upper_outliers(df_clean, "sex", "total_dur")
cat("\nРЕЗУЛЬТАТЫ УДАЛЕНИЯ ВЕРХНИХ ВЫБРОСОВ:\n")
##
## РЕЗУЛЬТАТЫ УДАЛЕНИЯ ВЕРХНИХ ВЫБРОСОВ:
outlier_summary <- df_clean %>%
group_by(sex) %>%
summarise(
было = n(),
среднее_было = round(mean(total_dur)/60, 1),
медиана_было = round(median(total_dur)/60, 1),
max_было = round(max(total_dur)/3600, 1) # максимальное значение в часах
) %>%
left_join(
df_no_outliers %>%
group_by(sex) %>%
summarise(
стало = n(),
среднее_стало = round(mean(total_dur)/60, 1),
медиана_стало = round(median(total_dur)/60, 1),
max_стало = round(max(total_dur)/3600, 1)
),
by = "sex"
) %>%
mutate(
удалено = было - стало,
процент_удалено = round(удалено / было * 100, 1),
изменение_среднего = среднее_стало - среднее_было,
изменение_медианы = медиана_стало - медиана_было,
изменение_max = max_стало - max_было
)
print(outlier_summary)
## # A tibble: 2 × 14
## sex было среднее_было медиана_было max_было стало среднее_стало
## <fct> <int> <dbl> <dbl> <dbl> <int> <dbl>
## 1 Ж 6051 144. 49.5 419. 5477 60.2
## 2 М 6446 125 41.9 973. 5995 57.6
## # ℹ 7 more variables: медиана_стало <dbl>, max_стало <dbl>, удалено <int>,
## # процент_удалено <dbl>, изменение_среднего <dbl>, изменение_медианы <dbl>,
## # изменение_max <dbl>
Ориентироваться на абсолютные цифры при асимметрии групп — правильнее, чем на проценты. Я удалил примерно равное количество экстремальных значений, а значит, сравннение “типичных мужчин” с “типичными женщинами” будет на равных основаниях.
# Разделяем данные по полу
men <- df_no_outliers$total_dur[df_no_outliers$sex == "М"]
women <- df_no_outliers$total_dur[df_no_outliers$sex == "Ж"]
# t-тест (параметрический) - ОСНОВНОЙ, так как дисперсии однородны
t_test_result <- t.test(men, women, var.equal = TRUE)
cat("\nt-тест (с равными дисперсиями):\n")
##
## t-тест (с равными дисперсиями):
print(t_test_result)
##
## Two Sample t-test
##
## data: men and women
## t = -2.2118, df = 11470, p-value = 0.027
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -291.64313 -17.58773
## sample estimates:
## mean of x mean of y
## 3458.694 3613.310
# U-тест Манна-Уитни (непараметрический) - для сравнения
wilcox_result <- wilcox.test(men, women)
cat("\nU-тест Манна-Уитни (для сравнения):\n")
##
## U-тест Манна-Уитни (для сравнения):
print(wilcox_result)
##
## Wilcoxon rank sum test with continuity correction
##
## data: men and women
## W = 16077432, p-value = 0.05508
## alternative hypothesis: true location shift is not equal to 0
p_value_t <- t_test_result$p.value
p_value_w <- wilcox_result$p.value
df_summary <- df_no_outliers %>%
group_by(sex) %>%
summarise(
count = n(),
mean_duration = mean(total_dur),
mean_minutes = round(mean(total_dur)/60, 1),
median_duration = median(total_dur),
median_minutes = round(median(total_dur)/60, 1),
sd_duration = sd(total_dur),
min_minutes = round(min(total_dur)/60, 1),
max_minutes = round(max(total_dur)/60, 1)
)
print(df_summary)
## # A tibble: 2 × 9
## sex count mean_duration mean_minutes median_duration median_minutes
## <fct> <int> <dbl> <dbl> <dbl> <dbl>
## 1 Ж 5477 3613. 60.2 2126 35.4
## 2 М 5995 3459. 57.6 1940 32.3
## # ℹ 3 more variables: sd_duration <dbl>, min_minutes <dbl>, max_minutes <dbl>
# кто смотрит дольше
women_mean <- df_summary$mean_minutes[df_summary$sex == "Ж"]
men_mean <- df_summary$mean_minutes[df_summary$sex == "М"]
women_median <- df_summary$median_minutes[df_summary$sex == "Ж"]
men_median <- df_summary$median_minutes[df_summary$sex == "М"]
diff_mean <- abs(women_mean - men_mean)
diff_median <- abs(women_median - men_median)
cat("\n Разница между группами:\n")
##
## Разница между группами:
cat(" • По среднему: женщины смотрят на", round(diff_mean, 1), "минут БОЛЬШЕ мужчин\n")
## • По среднему: женщины смотрят на 2.6 минут БОЛЬШЕ мужчин
cat(" • По медиане: женщины смотрят на", diff_median, "минут БОЛЬШЕ мужчин\n")
## • По медиане: женщины смотрят на 3.1 минут БОЛЬШЕ мужчин
# 7.1 Боксплот
p1 <- ggplot(df_no_outliers, aes(x = sex, y = total_dur/60, fill = sex)) +
geom_boxplot(alpha = 0.7,
width = 0.5,
outlier.size = 2,
outlier.alpha = 0.7,
color = "gray30") +
scale_fill_manual(values = c("М" = "#A4CCE8", "Ж" = "#FFC0CB")) +
theme_minimal() +
theme(
legend.position = "none",
plot.title = element_text(face = "bold", size = 14),
plot.subtitle = element_text(color = "gray40", size = 10),
axis.title = element_text(size = 11),
axis.text = element_text(size = 10),
panel.grid.major.x = element_blank(),
panel.grid.minor = element_blank()
) +
labs(title = "Распределение продолжительности просмотра по полу",
subtitle = paste("t-тест: p =", round(t_test_result$p.value, 4),
"| U-тест: p =", round(wilcox_result$p.value, 4)),
x = "Пол",
y = "Продолжительность просмотра (минуты)")
print(p1)
# 7.2 Столбчатая диаграмма со средними
p2 <- ggplot(df_summary, aes(x = sex, y = mean_minutes, fill = sex)) +
geom_col(alpha = 0.8,
width = 0.6,
color = "gray30",
linewidth = 0.3) +
scale_fill_manual(values = c("М" = "#A4CCE8", "Ж" = "#FFC0CB")) +
geom_text(aes(label = paste0(round(mean_minutes, 1), " мин")),
vjust = -0.8,
size = 4,
fontface = "bold") +
geom_errorbar(aes(ymin = mean_minutes - sd_duration/60,
ymax = mean_minutes + sd_duration/60),
width = 0.2,
linewidth = 0.8,
color = "gray30") +
theme_minimal() +
theme(
legend.position = "none",
plot.title = element_text(face = "bold", size = 14),
plot.subtitle = element_text(color = "gray40", size = 10),
axis.title = element_text(size = 11),
axis.text = element_text(size = 10),
panel.grid.major.x = element_blank(),
panel.grid.minor = element_blank(),
panel.grid.major.y = element_line(color = "gray90", linewidth = 0.3)
) +
labs(title = "Средняя продолжительность просмотра по полу",
subtitle = "Планки погрешности — стандартное отклонение",
x = "Пол",
y = "Средняя продолжительность (минуты)") +
scale_y_continuous(limits = c(0, max(df_summary$mean_minutes + df_summary$sd_duration/60) * 1.1),
expand = expansion(mult = c(0, 0.05)))
print(p2)
# 7.3 Столбчатая диаграмма с медианами
p3 <- ggplot(df_summary, aes(x = sex, y = median_minutes, fill = sex)) +
geom_col(alpha = 0.8,
width = 0.6,
color = "gray30",
linewidth = 0.3) +
scale_fill_manual(values = c("М" = "#A4CCE8", "Ж" = "#FFC0CB")) +
geom_text(aes(label = paste0(round(median_minutes, 1), " мин")),
vjust = -0.8,
size = 4,
fontface = "bold") +
theme_minimal() +
theme(
legend.position = "none",
plot.title = element_text(face = "bold", size = 14),
plot.subtitle = element_text(color = "gray40", size = 10),
axis.title = element_text(size = 11),
axis.text = element_text(size = 10),
panel.grid.major.x = element_blank(),
panel.grid.minor = element_blank(),
panel.grid.major.y = element_line(color = "gray90", linewidth = 0.3)
) +
labs(title = "Медианная продолжительность просмотра по полу",
subtitle = "Медиана устойчива к выбросам",
x = "Пол",
y = "Медианная продолжительность (минуты)") +
scale_y_continuous(limits = c(0, max(df_summary$median_minutes) * 1.2),
expand = expansion(mult = c(0, 0.05)))
print(p3)
Вопрос 1: Есть ли разница в продолжительности просмотра между мужчинами и женщинами?
• По t-тесту: ДА, различия значимы (p = 0.027)
• По U-тесту: НЕТ, различия не значимы (p = 0.055)
Учитывая, что распределение скошено, U-тест Манна-Уитни более надежен
Скорее всего, различий НЕТ
Вопрос 2: Кто в среднем смотрит контент дольше?
• Женщины: 60.2 минут (среднее), 35.4 минут (медиана)
• Мужчины: 57.6 минут (среднее), 32.3 минут (медиана)
Женщины смотрят на 2.6 минут БОЛЬШЕ по среднему, но эта разница очень мала ( 2.6 минут из 57.6 минут - около 4.5 %)
Вопрос 3: Значимы ли результаты?
• значим по t-тесту, но не значим по U-тесту
• U-тест Манна-Уитни более надежен для скошенных данных
Различия, скорее всего, НЕ ЗНАЧИМЫ
ОБЩИЙ ВЫВОД:
Несмотря на значимый t-тест (p = 0.027), U-тест Манна-Уитни (более надежный для наших данных) показывает отсутствие значимых различий (p = 0.055).
Даже если различия и есть, они крайне малы в реальных величинах - около 2.6 минут ( 4.5 %). Поэтому мы не можем уверенно утверждать, что пол влияет на продолжительность просмотра, следовательно подтвердилась нулевая гипотеза.
Нулевая гипотеза (H₀): Нет корреляционной связи между рейтингом контента на Кинопоиске и продолжительностью просмотра.
Альтернативная гипотеза (H₁): Существует статистически значимая корреляционная связь между рейтингом контента на Кинопоиске и продолжительностью просмотра.
Для анализа связи между двумя количественными переменными я буду использовать корреляционный анализ, поскольку обе переменные количественные:
rating_kp - рейтинг на Кинопоиске (от 0 до 10) total_dur - продолжительность просмотра в секундах
Так как мы уже знаем, что продолжительность просмотра имеет скошенное распределение и выбросы, я буду использовать ранговую корреляцию Спирмена.
# Подготовка данных
df_clean <- df %>%
filter(!is.na(rating_kp) & !is.na(total_dur) & total_dur > 0)
cat("Всего наблюдений после очистки:", nrow(df_clean), "\n\n")
## Всего наблюдений после очистки: 12637
Проверка нормальности
# Q-Q plots
qqnorm(df_clean$rating_kp, main = "Q-Q plot рейтингов")
qqline(df_clean$rating_kp, col = "red")
qqnorm(df_clean$total_dur, main = "Q-Q plot продолжительности")
qqline(df_clean$total_dur, col = "red")
par(mfrow = c(1, 1))
# Удаляем пропущенные значения для расчетов
rating_clean <- df_clean$rating_kp[!is.na(df_clean$rating_kp)]
# Создаем график
p1 <- ggplot(df_clean, aes(x = rating_kp)) +
# Гистограмма (автоматически пропускает NA)
geom_histogram(aes(y = after_stat(density)),
bins = 30,
fill = "#2E8B57",
alpha = 0.5,
color = "white",
linewidth = 0.2,
na.rm = TRUE) + # игнорируем NA
# Фактическая плотность
geom_density(color = "#D55E00",
linewidth = 1.2,
alpha = 0.8,
na.rm = TRUE) + # игнорируем NA
# Теоретическая нормальная кривая
stat_function(fun = dnorm,
args = list(mean = mean(rating_clean, na.rm = TRUE),
sd = sd(rating_clean, na.rm = TRUE)),
color = "#0072B2",
linetype = "dashed",
linewidth = 1) +
# Вертикальная линия для среднего
geom_vline(xintercept = mean(rating_clean, na.rm = TRUE),
color = "#333333",
linetype = "dotted",
linewidth = 0.8) +
theme_minimal() +
theme(
plot.title = element_text(face = "bold", size = 14),
plot.subtitle = element_text(color = "gray40", size = 10),
axis.title = element_text(size = 11),
axis.text = element_text(size = 10),
panel.grid.minor = element_blank(),
panel.grid.major = element_line(color = "gray90", linewidth = 0.3)
) +
labs(title = "Распределение рейтингов Кинопоиска",
subtitle = paste0("Красный — фактическое распределение, синий пунктир — нормальное",
"\n(пропущено ", sum(is.na(df_clean$rating_kp)), " значений)"),
x = "Рейтинг Кинопоиска",
y = "Плотность вероятности") +
# Аннотация со статистикой
annotate("text",
x = max(rating_clean, na.rm = TRUE) * 0.8,
y = max(density(rating_clean, na.rm = TRUE)$y) * 0.9,
label = paste0("Среднее = ", round(mean(rating_clean, na.rm = TRUE), 2),
"\nСт. откл. = ", round(sd(rating_clean, na.rm = TRUE), 2),
"\nN = ", length(rating_clean)),
hjust = 0,
size = 3.5,
color = "gray30")
print(p1)
# Для продолжительности (логарифмированная шкала)
# Сначала очищаем данные от NA и нулевых/отрицательных значений (log10 требует положительных)
dur_clean <- df_clean$total_dur[!is.na(df_clean$total_dur) & df_clean$total_dur > 0]
log_dur <- log10(dur_clean)
p2 <- ggplot(df_clean[!is.na(df_clean$total_dur) & df_clean$total_dur > 0, ],
aes(x = log10(total_dur))) +
# Гистограмма
geom_histogram(aes(y = after_stat(density)),
bins = 40, # чуть меньше для лучшей читаемости
fill = "#4682B4", # steelblue в HEX
alpha = 0.5,
color = "white",
linewidth = 0.2) +
# Фактическая плотность
geom_density(color = "#D55E00", # терракотовый
linewidth = 1.2,
alpha = 0.8) +
# Теоретическая нормальная кривая
stat_function(fun = dnorm,
args = list(mean = mean(log_dur),
sd = sd(log_dur)),
color = "#0072B2", # темно-синий
linetype = "dashed",
linewidth = 1) +
# Вертикальная линия для среднего
geom_vline(xintercept = mean(log_dur),
color = "#333333",
linetype = "dotted",
linewidth = 0.8) +
theme_minimal() +
theme(
plot.title = element_text(face = "bold", size = 14),
plot.subtitle = element_text(color = "gray40", size = 10),
axis.title = element_text(size = 11),
axis.text = element_text(size = 10),
panel.grid.minor = element_blank(),
panel.grid.major = element_line(color = "gray90", linewidth = 0.3)
) +
labs(title = "Распределение продолжительности просмотра (log10)",
subtitle = paste0("Красный — фактическое распределение, синий пунктир — нормальное"),
x = "log10 (продолжительность в секундах)",
y = "Плотность вероятности") +
# Аннотация со статистикой
annotate("text",
x = max(log_dur) * 0.8,
y = max(density(log_dur)$y) * 0.9,
label = paste0("Среднее (log10) = ", round(mean(log_dur), 3),
"\nСт. откл. (log10) = ", round(sd(log_dur), 3),
"\nМедиана (log10) = ", round(median(log_dur), 3),
"\nN = ", length(log_dur)),
hjust = 0,
size = 3.2,
color = "gray30") +
# Добавляем подписи на оригинальной шкале для удобства
scale_x_continuous(breaks = seq(floor(min(log_dur)), ceiling(max(log_dur)), by = 0.5),
labels = function(x) round(10^x, 0))
print(p2)
cat("ЧИСЛЕННЫЕ ХАРАКТЕРИСТИКИ РАСПРЕДЕЛЕНИЙ:\n")
## ЧИСЛЕННЫЕ ХАРАКТЕРИСТИКИ РАСПРЕДЕЛЕНИЙ:
# Для рейтинга
skew_rating <- skewness(df_clean$rating_kp, na.rm = TRUE)
kurt_rating <- kurtosis(df_clean$rating_kp, na.rm = TRUE)
cat("Рейтинг Кинопоиска:\n")
## Рейтинг Кинопоиска:
cat(" • Асимметрия (skewness):", round(skew_rating, 3))
## • Асимметрия (skewness): -1.887
if(abs(skew_rating) < 0.5) {
cat(" - близко к симметричному\n")
} else if(skew_rating > 0) {
cat(" - правосторонняя асимметрия\n")
} else {
cat(" - левосторонняя асимметрия\n")
}
## - левосторонняя асимметрия
cat(" • Эксцесс (kurtosis):", round(kurt_rating, 3))
## • Эксцесс (kurtosis): 10.465
if(abs(kurt_rating - 3) < 0.5) {
cat(" - близко к нормальному\n")
} else if(kurt_rating > 3) {
cat(" - островершинное (тяжелые хвосты)\n")
} else {
cat(" - плосковершинное\n")
}
## - островершинное (тяжелые хвосты)
# Для продолжительности
skew_dur <- skewness(df_clean$total_dur, na.rm = TRUE)
kurt_dur <- kurtosis(df_clean$total_dur, na.rm = TRUE)
cat("\nПродолжительность просмотра:\n")
##
## Продолжительность просмотра:
cat(" • Асимметрия (skewness):", round(skew_dur, 3))
## • Асимметрия (skewness): 51.889
if(abs(skew_dur) < 0.5) {
cat(" - близко к симметричному\n")
} else if(skew_dur > 0) {
cat(" - сильная правосторонняя асимметрия\n")
} else {
cat(" - левосторонняя асимметрия\n")
}
## - сильная правосторонняя асимметрия
cat(" • Эксцесс (kurtosis):", round(kurt_dur, 3))
## • Эксцесс (kurtosis): 3896.846
if(abs(kurt_dur - 3) < 0.5) {
cat(" - близко к нормальному\n")
} else if(kurt_dur > 3) {
cat(" - островершинное (очень тяжелые хвосты)\n")
} else {
cat(" - плосковершинное\n")
}
## - островершинное (очень тяжелые хвосты)
• Рейтинг Кинопоиска имеет распределение, близкое к нормальному
• Продолжительность просмотра имеет сильную правостороннюю асимметрию
• Есть явные выбросы в продолжительности просмотра
Исходя из этого, Для анализа лучше использовать корреляцию Спирмена
spearman_result <- cor.test(df_no_outliers$rating_kp,
df_no_outliers$total_dur,
method = "spearman",
exact = FALSE)
cat("\nКорреляция Спирмена:\n")
##
## Корреляция Спирмена:
cat(" ρ (rho) =", round(spearman_result$estimate, 4), "\n")
## ρ (rho) = -0.0012
cat(" S =", round(spearman_result$statistic, 0), "\n")
## S = 128117232907
cat(" p-value =", format(spearman_result$p.value, scientific = TRUE), "\n")
## p-value = 9.123827e-01
# Интерпретация силы связи
correlation_strength <- function(r) {
r <- abs(r)
if(r < 0.1) return("очень слабая")
else if(r < 0.2) return("слабая")
else if(r < 0.3) return("умеренная")
else if(r < 0.5) return("заметная")
else if(r < 0.7) return("сильная")
else return("очень сильная")
}
strength <- correlation_strength(spearman_result$estimate)
direction <- ifelse(spearman_result$estimate > 0, "положительная", "отрицательная")
Коэффициент корреляции Спирмена ρ = -0.0012 говорит о том, что:
• Связь между рейтингом и продолжительностью ПРАКТИЧЕСКИ ОТСУТСТВУЕТ
• Знак минус указывает на очень слабую отрицательную тенденцию, но значение настолько близко к нулю, что направление не имеет значения
p-value = 0.912 (> 0.05) означает, что:
• Вероятность получить такие данные случайно составляет 90%
• Результаты СТАТИСТИЧЕСКИ НЕ ЗНАЧИМЫ
# Очищаем данные от NA
df_scatter <- df_no_outliers[!is.na(df_no_outliers$rating_kp) &
!is.na(df_no_outliers$total_dur) &
df_no_outliers$total_dur > 0, ]
# Рассчитываем статистику для подписей
rho <- -0.0012
p_value <- 0.912
corr_text <- ifelse(p_value < 0.05, "Статистически значимая", "Статистически незначимая")
corr_color <- ifelse(p_value < 0.05, "#D55E00", "#0072B2")
p1 <- ggplot(df_scatter, aes(x = rating_kp, y = log10(total_dur))) +
# Точки
geom_point(alpha = 0.3,
size = 0.8,
color = "#1E3A5F",
shape = 16) +
# Линия регрессии с доверительным интервалом
geom_smooth(method = "lm",
color = corr_color,
fill = corr_color,
se = TRUE,
alpha = 0.15,
linewidth = 1) +
# Горизонтальная линия среднего
geom_hline(yintercept = mean(log10(df_scatter$total_dur)),
linetype = "dashed",
color = "gray40",
alpha = 0.7,
linewidth = 0.7) +
# Вертикальная линия среднего рейтинга
geom_vline(xintercept = mean(df_scatter$rating_kp, na.rm = TRUE),
linetype = "dashed",
color = "gray40",
alpha = 0.7,
linewidth = 0.7) +
# Тема без рамки
theme_minimal() +
theme(
plot.title = element_text(face = "bold", size = 15, hjust = 0.5),
plot.subtitle = element_text(color = corr_color,
size = 11,
hjust = 0.5,
face = "bold"),
axis.title = element_text(size = 12),
axis.text = element_text(size = 10),
panel.grid.minor = element_blank(),
panel.grid.major = element_line(color = "gray90", linewidth = 0.3),
plot.background = element_rect(fill = "white", color = NA),
panel.background = element_rect(fill = "white", color = NA), # убрана рамка панели
plot.margin = margin(20, 20, 20, 20)
) +
labs(title = "Связь рейтинга и продолжительности просмотра",
subtitle = paste0("ρ = ", round(rho, 4), " (p = ", round(p_value, 3), ") — ", corr_text),
x = "Рейтинг на Кинопоиске",
y = "log₁₀(продолжительность в секундах)") +
# Аннотация со статистикой в левом верхнем углу
annotate("rect",
xmin = 1.5, xmax = 3.0, # смещено влево
ymin = max(log10(df_scatter$total_dur)) * 0.82,
ymax = max(log10(df_scatter$total_dur)) * 0.95,
fill = "white",
alpha = 0.9,
color = NA) + # убран цвет рамки у прямоугольник
# Информация о распределении в левом нижнем углу
annotate("text",
x = min(df_scatter$rating_kp, na.rm = TRUE) + 0.2,
y = min(log10(df_scatter$total_dur)) + 0.1,
label = paste0("Ср. рейтинг: ", round(mean(df_scatter$rating_kp, na.rm = TRUE), 2),
"\nСр. длит: ", round(10^mean(log10(df_scatter$total_dur)), 0), " с"),
color = "gray40",
size = 3.2,
hjust = 0,
vjust = 0)
print(p1)
## `geom_smooth()` using formula = 'y ~ x'
# 5.2 Сглаженная кривая (для нелинейных зависимостей)
# Очищаем данные от NA
df_scatter <- df_no_outliers[!is.na(df_no_outliers$rating_kp) &
!is.na(df_no_outliers$total_dur) &
df_no_outliers$total_dur > 0, ]
# Рассчитываем статистику
mean_rating <- mean(df_scatter$rating_kp, na.rm = TRUE)
mean_log_dur <- mean(log10(df_scatter$total_dur), na.rm = TRUE)
p2 <- ggplot(df_scatter, aes(x = rating_kp, y = log10(total_dur))) +
# Точки с градиентом плотности
geom_point(alpha = 0.25,
size = 0.7,
color = "#1E5F3A", # темно-зеленый
shape = 16) +
# LOESS сглаживание с доверительным интервалом
geom_smooth(method = "loess",
color = "#0072B2", # синий
fill = "#0072B2",
se = TRUE,
alpha = 0.15,
linewidth = 1,
span = 0.75) + # параметр сглаживания
# Горизонтальная линия среднего
geom_hline(yintercept = mean_log_dur,
linetype = "dashed",
color = "gray40",
alpha = 0.7,
linewidth = 0.7) +
# Вертикальная линия среднего рейтинга
geom_vline(xintercept = mean_rating,
linetype = "dashed",
color = "gray40",
alpha = 0.7,
linewidth = 0.7) +
# Тема без рамки
theme_minimal() +
theme(
plot.title = element_text(face = "bold", size = 15, hjust = 0.5),
plot.subtitle = element_text(color = "gray40",
size = 11,
hjust = 0.5,
face = "plain"),
axis.title = element_text(size = 12),
axis.text = element_text(size = 10),
panel.grid.minor = element_blank(),
panel.grid.major = element_line(color = "gray90", linewidth = 0.3),
plot.background = element_rect(fill = "white", color = NA),
panel.background = element_rect(fill = "white", color = NA),
plot.margin = margin(20, 20, 20, 20)
) +
labs(title = "Связь рейтинга и продолжительности просмотра",
subtitle = "LOESS сглаживание — гориз. линия подтверждает отсутствие зависимости",
x = "Рейтинг на Кинопоиске",
y = "log10 (продолжительность в секундах)") +
# Аннотация со статистикой в левом верхнем углу
annotate("rect",
xmin = 1.5, xmax = 3.0,
ymin = max(log10(df_scatter$total_dur)) * 0.82,
ymax = max(log10(df_scatter$total_dur)) * 0.95,
fill = "white",
alpha = 0.9,
color = NA) +
# Информация о распределении в левом нижнем углу
annotate("text",
x = min(df_scatter$rating_kp, na.rm = TRUE) + 0.2,
y = min(log10(df_scatter$total_dur)) + 0.1,
label = paste0("Ср. рейтинг: ", round(mean_rating, 2),
"\nСр. длит: ", round(10^mean_log_dur, 0), " с",
"\nДов. интервал: 95%"),
color = "gray40",
size = 3.2,
hjust = 0,
vjust = 0) +
# Добавляем небольшую заметку о методе
annotate("text",
x = max(df_scatter$rating_kp, na.rm = TRUE) - 0.5,
y = max(log10(df_scatter$total_dur)) - 0.1,
label = "Метод: LOESS\n(локальная регрессия)",
color = "gray50",
size = 3,
hjust = 1,
vjust = 1,
alpha = 0.7)
print(p2)
## `geom_smooth()` using formula = 'y ~ x'
# Группировка рейтинга (без изменений)
df_no_outliers$rating_group <- cut(df_no_outliers$rating_kp,
breaks = c(0, 5, 6, 7, 8, 10),
labels = c("<5", "5-6", "6-7", "7-8", ">8"))
# Расчет глобальной медианы для горизонтальной линии
global_median <- median(df_no_outliers$total_dur, na.rm = TRUE) / 60
# Создание красивого графика
p3 <- ggplot(df_no_outliers, aes(x = rating_group, y = total_dur/60, fill = rating_group)) +
# Основные элементы
geom_boxplot(
alpha = 0.85, # немного прозрачности
width = 0.7, # ширина ящиков
outlier.shape = 21, # форма выбросов (кружок с заливкой)
outlier.size = 2, # размер выбросов
outlier.fill = "gray70", # заливка выбросов
outlier.alpha = 0.5, # прозрачность выбросов
notch = FALSE, # без вырезов для медианы
varwidth = FALSE, # одинаковая ширина ящиков
show.legend = FALSE
) +
# Медианная линия (без изменений)
geom_hline(
yintercept = global_median,
linetype = "dashed",
color = "#D55E00", # теплый оранжево-красный
alpha = 0.7,
linewidth = 0.8
) +
# Кастомная цветовая палитра (градиент зеленого)
scale_fill_manual(values = c(
"<5" = "#edf8e9", # очень светлый
"5-6" = "#c7e9c0", # светло-зеленый
"6-7" = "#74c476", # средний зеленый
"7-8" = "#31a354", # темно-зеленый
">8" = "#006d2c" # очень темный
)) +
# Ограничение оси Y (без изменений)
coord_cartesian(ylim = c(0, 150)) +
# Темы и оформление
theme_minimal(base_size = 14) +
theme(
# Текст
plot.title = element_text(face = "bold", size = 14, hjust = 0.5, margin = margin(b = 5)), # Уменьшено с 18 до 14
plot.subtitle = element_text(size = 13, hjust = 0.5, color = "gray30", margin = margin(b = 20)),
plot.caption = element_text(size = 10, color = "gray50", hjust = 1, margin = margin(t = 10)),
# Оси
axis.title.x = element_text(size = 14, margin = margin(t = 10)),
axis.title.y = element_text(size = 14, margin = margin(r = 10)),
axis.text.x = element_text(size = 12, color = "gray20"),
axis.text.y = element_text(size = 12, color = "gray20"),
axis.line = element_line(color = "gray70", linewidth = 0.3),
# Сетка
panel.grid.major = element_line(color = "gray90", linewidth = 0.3),
panel.grid.minor = element_blank(),
panel.background = element_rect(fill = "white", color = NA),
plot.background = element_rect(fill = "white", color = NA),
# Легенда (скрыта)
legend.position = "none",
# Рамка
plot.margin = margin(15, 20, 15, 20)
) +
# Подписи
labs(
title = "Распределение продолжительности просмотра",
subtitle = "Медианы практически одинаковы во всех группах",
x = "Рейтинг Кинопоиска",
y = "Продолжительность просмотра (минуты)",
caption = "Горизонтальная линия — общая медиана"
) +
# Добавление средних значений точками
stat_summary(
fun = mean,
geom = "point",
shape = 18, # ромбик
size = 3,
color = "#D55E00",
fill = "#D55E00",
show.legend = FALSE
)
print(p3)
# В начале файла добавьте загрузку библиотек
# Загружаем необходимые пакеты (без tidyverse)
library(ggplot2)
library(viridis) # для scale_fill_viridis
library(scales) # для comma
# Если пакеты не установлены, раскомментируйте:
# install.packages("ggplot2")
# install.packages("viridis")
# install.packages("scales")
# Тепловая карта плотности
p4 <- ggplot(df_no_outliers, aes(x = rating_kp, y = log10(total_dur))) +
# Основной слой с тепловой картой
geom_bin2d(
bins = 40,
alpha = 0.95
) +
# Цветовая шкала viridis
scale_fill_viridis(
option = "plasma",
direction = -1,
trans = "log",
labels = comma, # scales::comma
name = "Количество\nфильмов",
guide = guide_colorbar(
barwidth = 1.2,
barheight = 12,
title.position = "top",
title.hjust = 0.5
)
) +
# Линия регрессии
geom_smooth(
method = "lm",
se = TRUE,
color = "white",
fill = "gray90",
alpha = 0.3,
linetype = "dashed",
linewidth = 0.8,
formula = y ~ x
) +
# Средние линии
geom_vline(
xintercept = mean(df_no_outliers$rating_kp, na.rm = TRUE),
linetype = "dotted",
color = "white",
alpha = 0.5,
linewidth = 0.6
) +
geom_hline(
yintercept = mean(log10(df_no_outliers$total_dur), na.rm = TRUE),
linetype = "dotted",
color = "white",
alpha = 0.5,
linewidth = 0.6
) +
# Подписи
labs(
title = "Тепловая карта распределения",
subtitle = "Вертикальная структура подтверждает отсутствие зависимости",
x = "Рейтинг на Кинопоиске",
y = "log10(продолжительность)",
caption = "Пунктир — средние | Белая линия — регрессия"
) +
# Тема с уменьшенными размерами текста
theme_minimal(base_size = 9) +
theme(
# Текст
plot.title = element_text(face = "bold", size = 11, hjust = 0.5, margin = margin(b = 3)),
plot.subtitle = element_text(size = 9, hjust = 0.5, color = "gray30", margin = margin(b = 10)),
plot.caption = element_text(size = 7, color = "gray50", hjust = 1, margin = margin(t = 5)),
# Оси
axis.title.x = element_text(size = 9, margin = margin(t = 5)),
axis.title.y = element_text(size = 9, margin = margin(r = 5)),
axis.text.x = element_text(size = 8, color = "gray20"),
axis.text.y = element_text(size = 8, color = "gray20"),
axis.line = element_line(color = "gray70", linewidth = 0.3),
# Сетка
panel.grid.major = element_line(color = "gray85", linewidth = 0.3),
panel.grid.minor = element_blank(),
panel.background = element_rect(fill = "gray98", color = NA),
plot.background = element_rect(fill = "white", color = NA),
# Легенда
legend.position = "right",
legend.title = element_text(size = 8, face = "bold"),
legend.text = element_text(size = 7),
legend.background = element_rect(fill = "white", color = "gray80", linewidth = 0.3),
legend.margin = margin(5, 5, 5, 5),
legend.key.size = unit(0.8, "cm"),
# Рамка
plot.margin = margin(10, 15, 10, 15)
) +
# Настройка шкал
scale_x_continuous(
breaks = seq(0, 10, by = 1),
limits = c(0, 10),
expand = c(0.01, 0.01)
) +
scale_y_continuous(
breaks = seq(1, 5, by = 0.5),
expand = c(0.01, 0.01)
)
# Вывод графика
suppressWarnings(print(p4))
Вопрос 1: Связаны ли рейтинг контента и продолжительность просмотра?
НЕТ, связи не обнаружено (ρ = -0.0012, p = 0.912)
Вопрос 2: Какова сила и направление связи?
• Коэффициент корреляции практически равен нулю
• Связь полностью отсутствует
• Направление не имеет значения при таком коэффициенте
Вопрос 3: Значимы ли результаты?
НЕТ, результаты статистически не значимы (p = 0.912 > 0.05)
ОБЩИЙ ВЫВОД:
Несмотря на экстремальные характеристики распределений (асимметрия 51.9, эксцесс 3897), корреляционный анализ дает однозначный результат: связи между рейтингом фильма на Кинопоиске и продолжительностью его просмотра НЕ СУЩЕСТВУЕТ.
Пользователи смотрят контент одинаково долго независимо от того, имеет ли фильм высокий или низкий рейтинг. Рейтинг не является фактором, влияющим на вовлеченность зрителей, что подтверждает нулевую гипотезу.
Нулевая гипотеза (H₀): Пол пользователя и тип просматриваемого контента являются независимыми переменными. Предпочтения в выборе между фильмом и сериалом не зависят от пола пользователя.
Альтернативная гипотеза (H₁)): Между полом пользователя и типом просматриваемого контента существует статистически значимая связь. Предпочтения в выборе контента различаются в зависимости от пола.
Для проверки связи между двумя категориальными переменными я буду использовать критерий хи-квадрат χ² Пирсона. Этот выбор обусловлен тем, что обе переменные (“sex” и “content_type”) являются категориальными, и я хочу сравнить наблюдаемые частоты в группах с ожидаемыми, если бы связи не было.
Допущения теста хи-квадрат:
Наблюдения независимы.
Категории переменных взаимно исключают друг друга.
Ожидаемые частоты. Для достоверности теста не более 20% ячеек в таблице сопряженности должны иметь ожидаемую частоту меньше 5. Это условие будет проверено ниже.
contingency_table <- table(df$sex, df$content_type)
print("Таблица сопряженности:")
## [1] "Таблица сопряженности:"
print(contingency_table)
##
## film series
## Ж 4472 1579
## М 5184 1262
# Хи-квадрат тест
chi_result <- chisq.test(contingency_table)
print(chi_result)
##
## Pearson's Chi-squared test with Yates' continuity correction
##
## data: contingency_table
## X-squared = 75.091, df = 1, p-value < 2.2e-16
print("Ожидаемые частоты:")
## [1] "Ожидаемые частоты:"
print(chi_result$expected)
##
## film series
## Ж 4675.399 1375.601
## М 4980.601 1465.399
Все условия выполнены: все ожидаемые частоты > 5 и наблюдения не зависимы.
# Проценты по полу
props <- prop.table(contingency_table, 1) * 100
print("Проценты внутри каждой половой группы:")
## [1] "Проценты внутри каждой половой группы:"
print(props)
##
## film series
## Ж 73.90514 26.09486
## М 80.42197 19.57803
# Общие итоги
total_by_sex <- rowSums(contingency_table)
total_by_content <- colSums(contingency_table)
cat("\nВсего мужчин:", total_by_sex[1])
##
## Всего мужчин: 6051
cat("\nВсего женщин:", total_by_sex[2])
##
## Всего женщин: 6446
cat("\nВсего фильмов:", total_by_content[1])
##
## Всего фильмов: 9656
cat("\nВсего сериалов:", total_by_content[2])
##
## Всего сериалов: 2841
# Кто что предпочитает
cat("\n\n--- РЕЗУЛЬТАТЫ ---")
##
##
## --- РЕЗУЛЬТАТЫ ---
cat("\nМужчины:")
##
## Мужчины:
cat("\n Фильмы:", round(props[1,1], 1), "%")
##
## Фильмы: 73.9 %
cat("\n Сериалы:", round(props[1,2], 1), "%")
##
## Сериалы: 26.1 %
cat("\n\nЖенщины:")
##
##
## Женщины:
cat("\n Фильмы:", round(props[2,1], 1), "%")
##
## Фильмы: 80.4 %
cat("\n Сериалы:", round(props[2,2], 1), "%")
##
## Сериалы: 19.6 %
# Вывод о значимости
cat("\n\np-value =", chi_result$p.value)
##
##
## p-value = 4.494527e-18
p-value < 0.05 а значит результаты СТАТИСТИЧЕСКИ ЗНАЧИМЫ и можно отвергнуть нулевую гипотезу о независимости.
# Данные для графика
plot_df <- as.data.frame(contingency_table)
names(plot_df) <- c("Пол", "Тип_контента", "Количество")
# Проценты
plot_df <- plot_df %>%
group_by(Пол) %>%
mutate(Процент = Количество / sum(Количество) * 100)
# График
ggplot(plot_df, aes(x = Пол, y = Процент, fill = Тип_контента)) +
geom_col(position = "dodge") +
geom_text(aes(label = paste0(round(Процент, 1), "%")),
position = position_dodge(0.9), vjust = -0.3) +
labs(title = "Предпочтения контента по полу",
y = "Процент в группе (%)") +
theme_minimal()
Статистическая значимость:
Проведенный анализ с использованием критерия хи-квадрат Пирсона выявил статистически значимую связь между полом пользователя и типом просматриваемого контента (p < 0.05). Это позволяет отвергнуть нулевую гипотезу о независимости этих двух факторов.
Количественные различия в предпочтениях:
Женщины демонстрируют более высокую приверженность к просмотру фильмов (80.4% от всех просмотров), что на 6.5 процентных пункта выше, чем у мужчин
У мужчин фильмы составляют 73.9% от общего времени просмотра
Мужчины смотрят сериалы значительно чаще женщин (26.1% против 19.6%)
Разница в предпочтении сериалов между полами составляет 6.5 процентных пункта в пользу мужчин
Основной вывод:
Несмотря на то, что фильмы остаются основным типом контента для всех пользователей, мужчины проявляют больший интерес к сериалам, в то время как женщины отдают более сильное предпочтение фильмам. Это различие является статистически значимым, что отвергает нулевую гипотезу и подтверждает альтернативную.
Нулевая гипотеза (H₀): Предпочтения пользователей в типе просматриваемого контента (доля сериалов/фильмов) не зависят от уровня дохода. Распределение типов контента одинаково во всех доходных группах.
Альтернативная гипотеза (H₁): Предпочтения пользователей в типе просматриваемого контента зависят от уровня дохода. Существуют статистически значимые различия в том, какой контент предпочитают разные доходные группы.
Для проверки связи между двумя категориальными переменными я буду использовать критерий хи-квадрат χ² Пирсона. Этот выбор обусловлен тем, что обе переменные являются категориальными, и я хочу сравнить наблюдаемые частоты в группах с ожидаемыми, если бы связи не было.
Допущения теста хи-квадрат:
Наблюдения независимы.
Категории переменных взаимно исключают друг друга.
Ожидаемые частоты. Для достоверности теста не более 20% ячеек в таблице сопряженности должны иметь ожидаемую частоту меньше 5. Это условие будет проверено ниже.
# Смотрим распределение доходных групп
cat("Распределение по доходным группам:\n")
## Распределение по доходным группам:
print(table(df$income))
##
## income_0_20 income_150_inf income_20_40 income_40_60 income_60_90
## 288 14 6870 4006 1090
## income_90_150
## 247
# Создаем таблицу сопряженности доход × тип контента
income_content_table <- table(df$income, df$content_type)
cat("\nТаблица сопряженности (доход × тип контента):\n")
##
## Таблица сопряженности (доход × тип контента):
print(income_content_table)
##
## film series
## income_0_20 195 93
## income_150_inf 8 6
## income_20_40 5378 1492
## income_40_60 3086 920
## income_60_90 823 267
## income_90_150 178 69
# Убираем возможные пропущенные значения
df_clean <- df[!is.na(df$income) & !is.na(df$content_type), ]
cat("\nПосле удаления пропусков:", nrow(df_clean), "наблюдений\n")
##
## После удаления пропусков: 12515 наблюдений
# Добавляем simulate.p.value = TRUE для устранения предупреждения
chi_income_test <- chisq.test(income_content_table, simulate.p.value = TRUE, B = 2000)
# Проверка ожидаемых частот
cat("Проверка ожидаемых частот:")
## Проверка ожидаемых частот:
print(round(chi_income_test$expected, 1))
##
## film series
## income_0_20 222.5 65.5
## income_150_inf 10.8 3.2
## income_20_40 5307.2 1562.8
## income_40_60 3094.7 911.3
## income_60_90 842.0 248.0
## income_90_150 190.8 56.2
# Проверка условия применимости
expected_min <- min(chi_income_test$expected)
expected_less_5 <- sum(chi_income_test$expected < 5)
cat("Проверка условий теста:\n")
## Проверка условий теста:
cat("Минимальная ожидаемая частота:", round(expected_min, 1), "\n")
## Минимальная ожидаемая частота: 3.2
cat("Количество ячеек с ожидаемой частотой < 5:", expected_less_5, "\n")
## Количество ячеек с ожидаемой частотой < 5: 1
Поскольку в категории заработка “больше 150 тыс” данных меньше 5, я соединю эту категорию с категорией “доход 90-150”, тогда ожидаемая частота будет подходить под условие теста.
# Объединяем income_150_inf с income_90_150
df_merged <- df %>%
mutate(income_merged = case_when(
income == "income_150_inf" ~ "income_90_inf",
income == "income_90_150" ~ "income_90_inf",
TRUE ~ income
))
# Новая таблица
merged_table <- table(df_merged$income_merged, df_merged$content_type)
cat("ТАБЛИЦА ПОСЛЕ ОБЪЕДИНЕНИЯ ГРУПП:\n")
## ТАБЛИЦА ПОСЛЕ ОБЪЕДИНЕНИЯ ГРУПП:
print(merged_table)
##
## film series
## income_0_20 195 93
## income_20_40 5378 1492
## income_40_60 3086 920
## income_60_90 823 267
## income_90_inf 186 75
# Проверка ожидаемых частот после объединения
expected_merged <- chisq.test(merged_table)$expected
cat("\nОЖИДАЕМЫЕ ЧАСТОТЫ ПОСЛЕ ОБЪЕДИНЕНИЯ:\n")
##
## ОЖИДАЕМЫЕ ЧАСТОТЫ ПОСЛЕ ОБЪЕДИНЕНИЯ:
print(round(expected_merged, 1))
##
## film series
## income_0_20 222.5 65.5
## income_20_40 5307.2 1562.8
## income_40_60 3094.7 911.3
## income_60_90 842.0 248.0
## income_90_inf 201.6 59.4
# Хи-квадрат тест
chi_result <- chisq.test(merged_table)
cat("\nРЕЗУЛЬТАТЫ ТЕСТА ХИ-КВАДРАТ:\n")
##
## РЕЗУЛЬТАТЫ ТЕСТА ХИ-КВАДРАТ:
cat("X-squared =", round(chi_result$statistic, 3), "\n")
## X-squared = 26.403
cat("df =", chi_result$parameter, "\n")
## df = 4
cat("p-value =", format(chi_result$p.value, scientific = FALSE), "\n")
## p-value = 0.00002623663
p-value < 0.001, а значит результаты СТАТИСТИЧЕСКИ ЗНАЧИМЫ и можно отвергнуть нулевую гипотезу о независимости.
# Простой и надежный подход
# 1. Создаем dataframe из исходных данных
plot_df <- as.data.frame(merged_table)
names(plot_df) <- c("income", "content", "count")
# 2. Вручную создаем группы дохода
plot_df$income_group <- NA
for(i in 1:nrow(plot_df)) {
if(plot_df$income[i] == "income_0_20") plot_df$income_group[i] <- "0-20 тыс."
if(plot_df$income[i] == "income_20_40") plot_df$income_group[i] <- "20-40 тыс."
if(plot_df$income[i] == "income_40_60") plot_df$income_group[i] <- "40-60 тыс."
if(plot_df$income[i] == "income_60_90") plot_df$income_group[i] <- "60-90 тыс."
if(plot_df$income[i] %in% c("income_90_150", "income_150_inf")) plot_df$income_group[i] <- "более 90 тыс."
}
# 3. Проверяем, что все группы созданы
print("Проверка создания групп:")
## [1] "Проверка создания групп:"
print(table(plot_df$income_group))
##
## 0-20 тыс. 20-40 тыс. 40-60 тыс. 60-90 тыс.
## 2 2 2 2
# 4. Создаем новый dataframe с объединенными группами
# Сначала группируем по income_group и content, суммируем count
plot_df_agg <- aggregate(count ~ income_group + content, data = plot_df, FUN = sum)
# 5. Считаем проценты для каждой группы дохода
# Создаем отдельный dataframe с totals
totals <- aggregate(count ~ income_group, data = plot_df_agg, FUN = sum)
names(totals)[2] <- "total"
# Объединяем и считаем проценты
plot_df_agg <- merge(plot_df_agg, totals, by = "income_group")
plot_df_agg$percent <- plot_df_agg$count / plot_df_agg$total * 100
# 6. Устанавливаем порядок для графика
income_order <- c("0-20 тыс.", "20-40 тыс.", "40-60 тыс.",
"60-90 тыс.", "более 90 тыс.")
plot_df_agg$income_group <- factor(plot_df_agg$income_group, levels = income_order)
# 7. Проверяем финальные данные
print("Финальные данные для графика:")
## [1] "Финальные данные для графика:"
print(plot_df_agg)
## income_group content count total percent
## 1 0-20 тыс. film 195 288 67.70833
## 2 0-20 тыс. series 93 288 32.29167
## 3 20-40 тыс. film 5378 6870 78.28239
## 4 20-40 тыс. series 1492 6870 21.71761
## 5 40-60 тыс. film 3086 4006 77.03445
## 6 40-60 тыс. series 920 4006 22.96555
## 7 60-90 тыс. film 823 1090 75.50459
## 8 60-90 тыс. series 267 1090 24.49541
# 8. Строим график
p1 <- ggplot(plot_df_agg, aes(x = income_group, y = percent, fill = content)) +
geom_col(position = "dodge", width = 0.7, alpha = 0.9) +
scale_fill_manual(values = c("film" = "#3498db", "series" = "#e74c3c"),
labels = c("Фильмы", "Сериалы")) +
scale_y_continuous(limits = c(0, 100),
breaks = seq(0, 100, by = 20),
expand = expansion(mult = c(0, 0.05))) +
labs(title = "Предпочтения типа контента по доходным группам",
x = "Доходная группа",
y = "Доля внутри группы (%)",
fill = "Тип контента") +
theme_minimal() +
theme(
axis.text.x = element_text(angle = 45, hjust = 1, size = 11, face = "bold"),
legend.position = "bottom"
)
print(p1)
# Создаем данные для тепловой карты на основе merged_table
# 1. Берем исходные данные
heatmap_df <- as.data.frame(merged_table)
names(heatmap_df) <- c("income", "content", "count")
# 2. Создаем income_label с объединением групп
heatmap_df <- heatmap_df %>%
mutate(
income_label = case_when(
income == "income_0_20" ~ "0-20 тыс.",
income == "income_20_40" ~ "20-40 тыс.",
income == "income_40_60" ~ "40-60 тыс.",
income == "income_60_90" ~ "60-90 тыс.",
income %in% c("income_90_150", "income_150_inf") ~ "более 90 тыс.",
TRUE ~ as.character(income)
)
)
# 3. Группируем и считаем проценты для тепловой карты
heatmap_df <- heatmap_df %>%
group_by(income_label, content) %>%
summarise(count = sum(count), .groups = 'drop') %>%
group_by(income_label) %>%
mutate(
total = sum(count),
percent = count / total * 100
) %>%
ungroup()
# 4. Устанавливаем порядок для income_label
income_order <- c("0-20 тыс.", "20-40 тыс.", "40-60 тыс.",
"60-90 тыс.", "более 90 тыс.")
heatmap_df$income_label <- factor(heatmap_df$income_label, levels = income_order)
# 5. Проверяем данные
print("Данные для тепловой карты:")
## [1] "Данные для тепловой карты:"
print(heatmap_df)
## # A tibble: 10 × 5
## income_label content count total percent
## <fct> <fct> <int> <int> <dbl>
## 1 0-20 тыс. film 195 288 67.7
## 2 0-20 тыс. series 93 288 32.3
## 3 20-40 тыс. film 5378 6870 78.3
## 4 20-40 тыс. series 1492 6870 21.7
## 5 40-60 тыс. film 3086 4006 77.0
## 6 40-60 тыс. series 920 4006 23.0
## 7 60-90 тыс. film 823 1090 75.5
## 8 60-90 тыс. series 267 1090 24.5
## 9 <NA> film 186 261 71.3
## 10 <NA> series 75 261 28.7
# 6. Создаем тепловую карту
p3 <- ggplot(heatmap_df, aes(x = content, y = income_label, fill = percent)) +
geom_tile(color = "white", size = 1) +
geom_text(aes(label = paste0(round(percent, 1), "%")),
color = "white", size = 5, fontface = "bold") +
scale_fill_gradient(low = "#3498db", high = "#2c3e50",
name = "Доля (%)",
guide = guide_colorbar(barwidth = 15, barheight = 0.5)) +
scale_x_discrete(labels = c("film" = "Фильмы", "series" = "Сериалы")) +
labs(title = "Тепловая карта предпочтений",
subtitle = "Чем темнее цвет, тем выше доля просмотров",
x = "Тип контента",
y = "Доходная группа") +
theme_minimal() +
theme(
axis.text.x = element_text(size = 12, face = "bold", color = "#2c3e50"),
axis.text.y = element_text(size = 11, face = "bold"),
axis.title = element_text(size = 12, face = "bold"),
plot.title = element_text(size = 16, face = "bold", color = "#2c3e50"),
plot.subtitle = element_text(size = 12, color = "#7f8c8d"),
legend.position = "bottom",
legend.title = element_text(size = 11, face = "bold"),
legend.text = element_text(size = 10),
panel.grid = element_blank()
)
print(p3)
Предпочтения пользователей статистически значимо различаются в зависимости от уровня дохода, что отвергает нулевую гипотезу.
Результаты теста хи-квадрат:
X² = 26.403
df = 4
p-value = 0.0000262 < 0.001 → результаты высоко значимы. Вероятность случайного возникновения наблюдаемых различий составляет менее 0.003%.
props <- prop.table(merged_table, margin = 1) * 100
results_df <- data.frame(
Доходная_группа = rownames(props),
Фильмы = round(props[, 1], 1), # Берем первую колонку (фильмы)
Сериалы = round(props[, 2], 1), # Берем вторую колонку (сериалы)
Всего_просмотров = rowSums(merged_table),
stringsAsFactors = FALSE
)
names(results_df) <- c("Доходная группа", "Фильмы (%)", "Сериалы (%)", "Всего просмотров")
results_df$`Доходная группа` <- c(
"income_0_20" = "0-20 тыс. ₽",
"income_20_40" = "20-40 тыс. ₽",
"income_40_60" = "40-60 тыс. ₽",
"income_60_90" = "60-90 тыс. ₽",
"income_90_inf" = "более 90 тыс. ₽"
)[results_df$`Доходная группа`]
income_order <- c("0-20 тыс. ₽", "20-40 тыс. ₽", "40-60 тыс. ₽",
"60-90 тыс. ₽", "более 90 тыс. ₽")
results_df <- results_df[match(income_order, results_df$`Доходная группа`), ]
knitr::kable(results_df,
align = c('l', 'c', 'c', 'c'))
| Доходная группа | Фильмы (%) | Сериалы (%) | Всего просмотров | |
|---|---|---|---|---|
| income_0_20 | 0-20 тыс. ₽ | 67.7 | 32.3 | 288 |
| income_20_40 | 20-40 тыс. ₽ | 78.3 | 21.7 | 6870 |
| income_40_60 | 40-60 тыс. ₽ | 77.0 | 23.0 | 4006 |
| income_60_90 | 60-90 тыс. ₽ | 75.5 | 24.5 | 1090 |
| income_90_inf | более 90 тыс. ₽ | 71.3 | 28.7 | 261 |
📊 U-образная зависимость:
Группы с минимальным доходом (0-20 тыс.) и максимальным доходом (более 90 тыс.) демонстрируют повышенный интерес к сериалам (32.29% и 28.74% соответственно).
Группы со средним доходом (20-90 тыс.) показывают более сдержанный интерес к сериалам (21.72% - 24.50%).
Разрыв между группами:
Максимальная доля сериалов: 32.29% (0-20 тыс. ₽)
Минимальная доля сериалов: 21.72% (20-40 тыс. ₽)
Разница: 10.57 процентных пункта
Для продвижения сериалов наиболее перспективны две аудитории:
Пользователи с низким доходом (0-20 тыс. ₽) — 32.3% смотрят сериалы
Пользователи с высоким доходом (более 90 тыс. ₽) — 28.7% смотрят сериалы
Для продвижения фильмов лучше ориентироваться на средний класс:
Доход 20-40 тыс. ₽ — 78.3% смотрят фильмы
Доход 40-60 тыс. ₽ — 77.0% смотрят фильмы
U-образная зависимость говорит о том, что стратегия продвижения контента должна быть разной для разных доходных групп:
Бедные и богатые — делаем упор на сериалы
Средний класс — делаем упор на фильмы