Разведочный анализ в R
2024-03-27
Провести комплексную проверку на нормальность переменной Sepal.Length из набора `iris.
library(DataExplorer)i<-irisintroduce(iris)## rows columns discrete_columns continuous_columns all_missing_columns
## 1 150 5 1 4 0
## total_missing_values complete_rows total_observations memory_usage
## 1 0 150 750 7976plot_intro(iris)
summary(iris$Sepal.Length)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 4.300 5.100 5.800 5.843 6.400 7.900library(psych)psych::describe(iris$Sepal.Length)## vars n mean sd median trimmed mad min max range skew kurtosis se
## X1 1 150 5.84 0.83 5.8 5.81 1.04 4.3 7.9 3.6 0.31 -0.61 0.07psych::describe(iris$Sepal.Length)## vars n mean sd median trimmed mad min max range skew kurtosis se
## X1 1 150 5.84 0.83 5.8 5.81 1.04 4.3 7.9 3.6 0.31 -0.61 0.07library(ggplot2)##
## Attaching package: 'ggplot2'## The following objects are masked from 'package:psych':
##
## %+%, alphalibrary(DataExplorer)plot_intro(psych::describe(iris$Sepal.Length))
hist(iris$Sepal.Length, col="steelblue", xlab = "длина", main="Чашелистник")
library(ggpubr)ggpubr::ggqqplot(iris, x = "Sepal.Length", color = "steelblue")
boxplot(iris$Sepal.Length, col="steelblue")
iris$z_Sepal.Length <- scale(iris$Sepal.Length)
hist(iris$z_Sepal.Length)
summary(iris$z_Sepal.Length)## V1
## Min. :-1.86378
## 1st Qu.:-0.89767
## Median :-0.05233
## Mean : 0.00000
## 3rd Qu.: 0.67225
## Max. : 2.48370which(iris$z_Sepal.Length > 2.0)## [1] 106 118 119 123 132 136iris[c(106,118,119,123,132,136),]$Sepal.Length## [1] 7.6 7.7 7.7 7.7 7.9 7.7library(outliers)##
## Attaching package: 'outliers'## The following object is masked from 'package:psych':
##
## outliertest <- grubbs.test(iris$Sepal.Length, type=10)
test##
## Grubbs test for one outlier
##
## data: iris$Sepal.Length
## G = 2.48370, U = 0.95832, p-value = 0.916
## alternative hypothesis: highest value 7.9 is an outlierlibrary(ggplot2)
ggplot(iris, aes(x=Sepal.Length)) +
geom_density(fill="lightblue")+
theme_bw()
library(DescTools)##
## Attaching package: 'DescTools'## The following objects are masked from 'package:psych':
##
## AUC, ICC, SDDesc(iris$Sepal.Length, plotit = TRUE)## ------------------------------------------------------------------------------
## iris$Sepal.Length (numeric)
##
## length n NAs unique 0s mean meanCI'
## 150 150 0 35 0 5.843 5.710
## 100.0% 0.0% 0.0% 5.977
##
## .05 .10 .25 median .75 .90 .95
## 4.600 4.800 5.100 5.800 6.400 6.900 7.255
##
## range sd vcoef mad IQR skew kurt
## 3.600 0.828 0.142 1.038 1.300 0.309 -0.606
##
## lowest : 4.3, 4.4 (3), 4.5, 4.6 (4), 4.7 (2)
## highest: 7.3, 7.4, 7.6, 7.7 (4), 7.9
##
## ' 95%-CI (classic)
Skew(iris$Sepal.Length, na.rm=TRUE, conf.level = 0.95, ci.type = "classic")## skewness lwr.ci upr.ci
## 0.3086407 -0.3804196 0.3804196Kurt(iris$Sepal.Length, na.rm=TRUE, method=2, conf.level = 0.95, ci.type = "classic")## kurtosis lwr.ci upr.ci
## -0.5520640 -0.7714087 0.7714087shapiro.test(iris$Sepal.Length)##
## Shapiro-Wilk normality test
##
## data: iris$Sepal.Length
## W = 0.97609, p-value = 0.01018library(nortest)nortest::lillie.test(iris$Sepal.Length)##
## Lilliefors (Kolmogorov-Smirnov) normality test
##
## data: iris$Sepal.Length
## D = 0.088654, p-value = 0.005788Провести одномерный анализ по переменным V12, V13, V15. Сделать двумерный анализ по региону и возрасту.
library(haven)
База_КлимРиск_2023_ОБЩАЯ_от_31_10_23_с_правильными_типами_для_отчета <- read_sav("База_КлимРиск_2023_ОБЩАЯ от 31.10.23_с правильными типами для отчета.sav")df<-База_КлимРиск_2023_ОБЩАЯ_от_31_10_23_с_правильными_типами_для_отчетаlibrary(psych)
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, uniondf <- haven::as_factor(df)
table(df$V12)##
## Выше Ниже Не изменилась
## 190 437 189
## Другое Затрудняюсь ответить
## 0 0prop.table(table(df$V12))*100##
## Выше Ниже Не изменилась
## 23.28431 53.55392 23.16176
## Другое Затрудняюсь ответить
## 0.00000 0.00000df <- haven::as_factor(df)
table(df$V13)##
## Стала холоднее Стала теплее Не изменилась
## 575 149 166prop.table(table(df$V13))*100##
## Стала холоднее Стала теплее Не изменилась
## 64.60674 16.74157 18.65169df <- haven::as_factor(df)
table(df$V15)##
## Стала холоднее Стала теплее Не изменилась
## 586 153 154prop.table(table(df$V15))*100##
## Стала холоднее Стала теплее Не изменилась
## 65.62150 17.13326 17.24524##Сделать двумерный анализ по региону и возрасту.
library(gtsummary)table<-df %>%
select(Region, age_cats3) table %>%
tbl_summary(missing="no")| Characteristic | N = 9131 |
|---|---|
| Регион | |
| Алтайский край | 274 (30%) |
| Республика Алтай | 446 (49%) |
| Республика Тыва | 192 (21%) |
| Монголия | 0 (0%) |
| Возраст | |
| До 30 лет | 160 (18%) |
| 31-49 лет | 424 (48%) |
| 50 лет и старше | 303 (34%) |
| 1 n (%) | |
table %>%
tbl_summary(by=Region, missing="no") %>%
modify_header(
update = list(
label ~ "Характеристика"
)
)## 1 observations missing `Region` have been removed. To include these observations, use `forcats::fct_na_value_to_level()` on `Region` column before passing to `tbl_summary()`.| Характеристика | Алтайский край, N = 2741 | Республика Алтай, N = 4461 | Республика Тыва, N = 1921 | Монголия, N = 01 |
|---|---|---|---|---|
| Возраст | ||||
| До 30 лет | 35 (13%) | 75 (17%) | 50 (27%) | 0 (NA%) |
| 31-49 лет | 114 (43%) | 215 (49%) | 94 (51%) | 0 (NA%) |
| 50 лет и старше | 115 (44%) | 147 (34%) | 41 (22%) | 0 (NA%) |
| 1 n (%) | ||||
table %>%
tbl_cross(row = Region, col = age_cats3, percent = "cell", missing="no", ) %>%
add_p(test="chisq.test") %>%
modify_header(
update = list(
label ~ "**Характеристика**" #
)
) %>%
modify_caption("**Взаимосвязь регина и возраста** (N = {N})") ## FALSE observations with missing data have been removed.| Характеристика | Возраст | Total | p-value1 | ||
|---|---|---|---|---|---|
| До 30 лет | 31-49 лет | 50 лет и старше | |||
| Регион | <0.001 | ||||
| Алтайский край | 35 (4.0%) | 114 (13%) | 115 (13%) | 264 (30%) | |
| Республика Алтай | 75 (8.5%) | 215 (24%) | 147 (17%) | 437 (49%) | |
| Республика Тыва | 50 (5.6%) | 94 (11%) | 41 (4.6%) | 185 (21%) | |
| Монголия | 0 (0%) | 0 (0%) | 0 (0%) | 0 (0%) | |
| Total | 160 (18%) | 423 (48%) | 303 (34%) | 886 (100%) | |
| 1 Pearson’s Chi-squared test | |||||
age<- as.data.frame(sjmisc::flat_table(df, age_cats3, margin = "cell", digits = 1))
age_Region<-as.data.frame(sjmisc::flat_table(df, age_cats3, Region, margin = "row", digits = 1))library(ggplot2)
library(scales)##
## Attaching package: 'scales'## The following objects are masked from 'package:psych':
##
## alpha, rescaleggplot(data=age, aes(x=age_cats3, y=Freq)) +
geom_bar(stat="identity", fill="steelblue", width = 0.6)+
geom_text(aes(label=number(Freq,accuracy=0.1, decimal.mark = ",")), family = "Ink Free", vjust=0.5, hjust=-0.5, color="steelblue", size=4)+
coord_flip()+
theme_void()+
expand_limits(y=c(0,50))+
theme(axis.text.y = element_text(size = 10, family = "Ink Free"))+
scale_x_discrete(limits=rev)
ggplot(data=age_Region, aes(x=age_cats3, y=Freq, fill=Region)) +
geom_bar(stat="identity", position="dodge", width = 0.6)+
geom_text(aes(label=number(Freq,accuracy=0.1, decimal.mark = ",")), position = position_dodge(width = 0.6), hjust=-0.3, family = "Ink Free", size=4)+
coord_flip()+
theme_void()+
expand_limits(y=c(0,100))+
theme(axis.text.y = element_text(size = 10, family = "Ink Free"))+
scale_x_discrete(limits=rev)+
theme(legend.position="bottom")+
theme(legend.title=element_blank())
ggstatsplot::ggpiestats(df, Region, age_cats3)
chisq.test(df$Region, df$age_cats3)##
## Pearson's Chi-squared test
##
## data: df$Region and df$age_cats3
## X-squared = 28.156, df = 4, p-value = 1.16e-05Проанализировать переменные с множественным выбором V14 и V16, также сделать двумерный анализ по региону.
library(questionr)##
## Attaching package: 'questionr'## The following object is masked from 'package:psych':
##
## describelibrary(dplyr)
V14<-df %>%
select(contains("V14")) %>%
select_if(is.factor)library(questionr)
library(dplyr)
V14tab<-V14 %>%
multi.table(true.codes=list("да"), freq=TRUE) %>%
as.data.frame() %>%
rename(Pct=`%multi`) %>%
arrange(desc(Pct))
V14tab## n Pct
## V14_2 436 47.8
## V14_1 428 46.9
## V14_10 381 41.7
## V14_5 310 34.0
## V14_4 266 29.1
## V14_7 211 23.1
## V14_3 102 11.2
## V14_8 64 7.0
## V14_6 58 6.4
## V14_9 42 4.6
## V14_88 23 2.5library(sjlabelled)##
## Attaching package: 'sjlabelled'## The following object is masked from 'package:dplyr':
##
## as_label## The following objects are masked from 'package:haven':
##
## as_factor, read_sas, read_spss, read_stata, write_sas, zap_labels## The following object is masked from 'package:ggplot2':
##
## as_labelV14labs<-df %>%
select(contains("V14")) %>%
select_if(is.factor) %>%
get_label()
V14labs## V14_1
## "Резкие перепады температуры (от холода к теплу и наоборот)"
## V14_2
## "Затяжные морозы, увеличение периодов аномального холода"
## V14_3
## "Частые оттепели"
## V14_4
## "Увеличение гололеда на дорогах"
## V14_5
## "Обильные снегопады, увеличение снежного покрова"
## V14_6
## "Сход снежных лавин"
## V14_7
## "Увеличение количества пасмурных дней, нехватка солнца"
## V14_8
## "Раннее таяние и более позднее образование речного льда"
## V14_9
## "Участились ледовые зажоры, наледи на реках"
## V14_10
## "Усилились ветра, метели и снежные наносы"
## V14_88
## "Другое"V14tab$n<-V14labsV14tab## n Pct
## V14_2 Резкие перепады температуры (от холода к теплу и наоборот) 47.8
## V14_1 Затяжные морозы, увеличение периодов аномального холода 46.9
## V14_10 Частые оттепели 41.7
## V14_5 Увеличение гололеда на дорогах 34.0
## V14_4 Обильные снегопады, увеличение снежного покрова 29.1
## V14_7 Сход снежных лавин 23.1
## V14_3 Увеличение количества пасмурных дней, нехватка солнца 11.2
## V14_8 Раннее таяние и более позднее образование речного льда 7.0
## V14_6 Участились ледовые зажоры, наледи на реках 6.4
## V14_9 Усилились ветра, метели и снежные наносы 4.6
## V14_88 Другое 2.5V14Regiontab<-cross.multi.table(V14, df$Region, true.codes=list("да"), freq=TRUE)
V14Regiontab<-as_tibble(V14Regiontab)
V14Regiontab$n<-V14labs
V14Regiontab<-V14Regiontab %>%
relocate(n)V14Regiontab## # A tibble: 11 × 4
## n `Алтайский край` `Республика Алтай` `Республика Тыва`
## <chr> <dbl> <dbl> <dbl>
## 1 Резкие перепады темпер… 63.5 43.3 31.8
## 2 Затяжные морозы, увели… 25.5 63 44.3
## 3 Частые оттепели 17.9 9.4 5.7
## 4 Увеличение гололеда на… 23.4 29.6 36.5
## 5 Обильные снегопады, ув… 22.6 39.7 37
## 6 Сход снежных лавин 2.2 9 6.2
## 7 Увеличение количества … 17.5 26.7 22.9
## 8 Раннее таяние и более … 7.7 8.1 3.6
## 9 Участились ледовые заж… 1.8 6.7 3.6
## 10 Усилились ветра, метел… 34.7 47.3 39.1
## 11 Другое 0 4.5 1.6install.packages("extrafont")## Installing package into '/cloud/lib/x86_64-pc-linux-gnu-library/4.3'
## (as 'lib' is unspecified)library(extrafont)## Registering fonts with Rfont_import()## Importing fonts may take a few minutes, depending on the number of fonts and the speed of the system.
## Continue? [y/n]## Exiting.library(extrafont)
font_import()## Importing fonts may take a few minutes, depending on the number of fonts and the speed of the system.
## Continue? [y/n]## Exiting.loadfonts(device = "win", quiet = TRUE)## Warning in loadfonts_win(quiet = quiet): OS is not Windows. No fonts registered
## with windowsFonts().library(dplyr)
library(gtsummary)
library(ggplot2)
library(scales)library(stringr)
V14tab %>%
ggplot(aes(x=n, y=Pct, fill=Pct)) + # создадим типовой график
geom_bar(stat="identity", width = 0.6)+ # добавим geom
geom_text(aes(label=number(Pct,accuracy=0.1, decimal.mark = ",")), hjust=-0.3, size=4)+ # поработаем с подписями данных
expand_limits(y=c(0,80))+# увеличим лимит оси
coord_flip()+ # перевернем график в горизонтальное положение
theme_void()+ # добавим минималистичную тему
theme(axis.text.y = element_text(size = 10, hjust=1))+ # установим шрифт для подписей
theme(legend.position="none")+ # уберем легенду
scale_fill_gradient(low = "blue", high = "green")+ # установим цвета градиента для столбцов
scale_x_discrete(labels = function(x) str_wrap(x, width = 30))
V16<-df %>%
select(contains("V16")) %>%
select_if(is.factor)library(questionr)
library(dplyr)
V16tab<-V16 %>%
multi.table(true.codes=list("да"), freq=TRUE) %>%
as.data.frame() %>%
rename(Pct=`%multi`) %>%
arrange(desc(Pct))
V16tab## n Pct
## V16_7 406 44.5
## V16_3 392 42.9
## V16_1 326 35.7
## V16_2 290 31.8
## V16_8 221 24.2
## V16_5 173 18.9
## V16_9 142 15.6
## V16_6 80 8.8
## V16_4 57 6.2
## V16_88 34 3.7library(sjlabelled)
V16labs<-df %>%
select(contains("V16")) %>%
select_if(is.factor) %>%
get_label()
V16labs## V16_1
## "Увеличение количества засушливых дней, без осадков"
## V16_2
## "Увеличение периодов аномальной жары"
## V16_3
## "Большое количество осадков, сильных дождей"
## V16_4
## "Увеличение количества камнепадов и оползней в горах"
## V16_5
## "Усилился паводок, разливы рек, затопление мест, которые раньше не затапливались"
## V16_6
## "Таяние вечной мерзлоты, выход грунтовых вод на поверхность"
## V16_7
## "Сильные ветры, штормы"
## V16_8
## "Рост количества насекомых, комаров, мошки"
## V16_9
## "Рост количества насекомых-вредителей, угрожающих сельскохозяйственным культурам, хвойным лесам"
## V16_88
## "Другое"V16tab$n<-V16labsV16tab## n
## V16_7 Увеличение количества засушливых дней, без осадков
## V16_3 Увеличение периодов аномальной жары
## V16_1 Большое количество осадков, сильных дождей
## V16_2 Увеличение количества камнепадов и оползней в горах
## V16_8 Усилился паводок, разливы рек, затопление мест, которые раньше не затапливались
## V16_5 Таяние вечной мерзлоты, выход грунтовых вод на поверхность
## V16_9 Сильные ветры, штормы
## V16_6 Рост количества насекомых, комаров, мошки
## V16_4 Рост количества насекомых-вредителей, угрожающих сельскохозяйственным культурам, хвойным лесам
## V16_88 Другое
## Pct
## V16_7 44.5
## V16_3 42.9
## V16_1 35.7
## V16_2 31.8
## V16_8 24.2
## V16_5 18.9
## V16_9 15.6
## V16_6 8.8
## V16_4 6.2
## V16_88 3.7V16Regiontab<-cross.multi.table(V16, df$Region, true.codes=list("да"), freq=TRUE)
V16Regiontab<-as_tibble(V16Regiontab)
V16Regiontab$n<-V16labs
V16Regiontab<-V16Regiontab %>%
relocate(n)V16Regiontab## # A tibble: 10 × 4
## n `Алтайский край` `Республика Алтай` `Республика Тыва`
## <chr> <dbl> <dbl> <dbl>
## 1 Увеличение количества … 23.7 45.3 30.7
## 2 Увеличение периодов ан… 39.1 31.4 22.4
## 3 Большое количество оса… 60.2 34.1 39.1
## 4 Увеличение количества … 0.7 8.3 9.4
## 5 Усилился паводок, разл… 6.9 27.6 16.1
## 6 Таяние вечной мерзлоты… 2.6 14.8 3.6
## 7 Сильные ветры, штормы 29.2 52.2 48.4
## 8 Рост количества насеко… 26.6 26.5 15.6
## 9 Рост количества насеко… 20.8 15.9 7.3
## 10 Другое 0.4 5.6 4.2library(stringr)
V16tab %>%
ggplot(aes(x=n, y=Pct, fill=Pct)) +
geom_bar(stat="identity", width = 1)+
geom_text(aes(label=number(Pct,accuracy=0.1, decimal.mark = ",")), hjust=-0.5, size=3)+
expand_limits(y=c(1))+
coord_flip()+
theme_void()+
theme(axis.text.y = element_text(size = 10, hjust=1))+
theme(legend.position="none")+
scale_fill_gradient(low = "blue", high = "green")+
scale_x_discrete(labels = function(x) str_wrap(x, width = 40))