微信公众号搜"智元新知"关注
微信扫一扫可直接关注哦!

ggplot:如何“纠正”图中的代表性峰值 但是,事实证明,有一天中的某一天恰好具有较少的数据点因此,实际数据集如下所示:我一直在做的事情依靠手段

如何解决ggplot:如何“纠正”图中的代表性峰值 但是,事实证明,有一天中的某一天恰好具有较少的数据点因此,实际数据集如下所示:我一直在做的事情依靠手段

我有沿日期时间(日期和小时:分钟:秒)的百分比得分数据。我想以图形方式“纠正” /突出显示不具有代表性的数据点。

背景

我有关于人们每天如何以0到0的连续等级来评价自己的幸福水平的数据,其中0表示“极度不高兴”,1表示“极度不高兴”。我问很多人,并希望随着时间的流逝变得“团体中的幸福”。

数据

library(tidyverse)
library(lubridate)

set.seed(1234)

original_df <- 
  seq(as.POSIXct('2020-09-01',tz = "UTC"),as.POSIXct('2020-09-15',by="1 mins") %>%
  sample(15000,replace = T) %>%
  as_tibble %>%
  rename(date_time = value) %>%
  mutate(date = date(date_time)) %>%
  add_column(score = runif(15000))

original_df

## # A tibble: 15,000 x 3
##  date_time           date       score
##    <dttm>              <date>     <dbl>
##  1 2020-09-06 04:11:00 2020-09-06 0.683
##  2 2020-09-06 13:35:00 2020-09-06 0.931
##  3 2020-09-05 23:21:00 2020-09-05 0.121
##  4 2020-09-06 14:45:00 2020-09-06 0.144
##  5 2020-09-07 09:15:00 2020-09-07 0.412
##  6 2020-09-01 10:22:00 2020-09-01 0.564
##  7 2020-09-11 14:00:00 2020-09-11 0.960
##  8 2020-09-08 13:24:00 2020-09-08 0.845
##  9 2020-09-01 15:33:00 2020-09-01 0.225
## 10 2020-09-09 19:27:00 2020-09-09 0.815
## # ... with 14,990 more rows

但是,事实证明,有一天中的某一天恰好具有较少的数据点。因此,实际数据集如下所示:

actual_df <- 
  original_df %>%
  filter(date %in% as_date("2020-09-10")) %>%
  group_by(date) %>%
  slice_sample(n = 15) %>%
  ungroup %>%
  bind_rows(original_df %>% filter(!date %in% as_date("2020-09-10")))

> actual_df %>% count(date)

## # A tibble: 14 x 2
##    date           n
##    <date>     <int>
##  1 2020-09-01  1073
##  2 2020-09-02  1079
##  3 2020-09-03  1118
##  4 2020-09-04  1036
##  5 2020-09-05  1025
##  6 2020-09-06  1089
##  7 2020-09-07  1040
##  8 2020-09-08  1186
##  9 2020-09-09  1098
## 10 2020-09-10    15 ## <- this day has less data 
## 11 2020-09-11  1095
## 12 2020-09-12  1051
## 13 2020-09-13  1037
## 14 2020-09-14  1034

绘制此数据随时间变化

我一直在做的事情依靠手段

我将每天视为一个因素,并获得每日平均值。从统计学上讲,此解决方案可能远非理想,如下面的@BrianLang所述。但是,现在这是我选择的方法

library(emmeans)

model_fit <- 
  actual_df %>%
  mutate(across(date,factor)) %>%
  lm(score ~ date,data = .)

emmeans_fit_data <- emmeans(model_fit,~ date,CIs = TRUE)

emmeans_fit_data %>%
  as_tibble %>%
  ggplot(data = .,aes(x = date,y = emmean)) +
  geom_line(color = "#1a476f",group = 1,lwd = 1) +
  geom_errorbar(aes(ymin = lower.CL,ymax = upper.CL),alpha = 0.5,color = "#90353b",width = 0.2) +
  geom_text(aes(label = paste0(round(100*emmean,1),"%"),color = "90353b"),vjust = -4,hjust = 0.5,size = 3.5) +
  geom_point(color = "1a476f") +
  scale_y_continuous(labels = function(x) paste0(100*x,"%")) +
  ylab("Level of Happiness") +
  xlab("Date") +
  ggtitle("Mood Over Time") +
  theme(plot.title = element_text(hjust = 0.5,size = 14),axis.text.x=element_text(angle = -60,hjust = 0),axis.title.x = element_blank(),legend.title = element_blank(),plot.caption = element_text(hjust = 0,size = 8),legend.position = "none")

enter image description here



但是后来我在2020-09-10达到了峰值,这仅仅是由于数据点数量少。 一种图形化的解决方案是做一些事情,例如将有问题的虚线划线并“完成”具有足够数据点的情况。也许基于前一天和后一天的平均值?我不想摆脱真实的数据,但是想以图形方式突出显示这是不具有代表性的,真实的价值应该更接近前后。我当时认为使用虚线是一种合理的图形解决方案。

dashed



否则,我希望可以使用ggplot的平滑方法来对“按时间”数据进行建模/绘图的方法不同,这将为我提供一条更平滑的趋势线和一条置信带解决有问题的一天。但是我知道这可能超出了这个问题的范围,因此我只是将其添加为补充说明;如果有人想提出一个基于不同模型的解决方案,而不是单纯的图形校正。但我都会为此而感激。

解决方法

不想进入时间序列模型,您可以想象使用受限三次样条曲线变换时间变量。

我需要更改一些代码,以便避免安装某些软件包的最新版本;-)。

请注意,我更改了一些变量名,因为date是一个函数名,不应同时用作变量名。

library(chron)

## added a numeric version of your date variable.
actual_df <- original_df %>%
 filter(datez %in% lubridate::date("2020-09-10")) %>%
 sample_n(size = 15) %>%
 group_by(datez) %>%
 ungroup %>%
 bind_rows(original_df %>% filter(!datez %in% lubridate::date("2020-09-10"))) %>%
 mutate(num_date = as.numeric(datez))
## How many knots across the dates do you want?
number_of_knots = 15

## This is to make sure that visreg is passed the actual knot locations! RMS::RCS does not store them in the model fits. 
knots <- paste0("c(",paste0(attr(rms::rcs(actual_df$num_date,number_of_knots),"parms"),collapse = ","),")") 

## We can construct the formula early.
formula <- as.formula(paste("score ~ rms::rcs(num_date,",knots,")"))

## fit the model as a gaussian glm and pass it to visreg for it's prediction function. This will give you predicted means and 95% CI for that mean. Then I convert the numeric dates back to real dates. 
glm_rcs <- glm(data = actual_df,formula = formula,family = "gaussian") %>% visreg::visreg(plot = F) %>% .$fit %>%
 mutate(date_date = chron::as.chron(num_date) %>% as.POSIXct())

## plot it!
ggplot(data = glm_rcs,aes(date_date,y = visregFit)) + 
 geom_ribbon(aes(ymin = visregLwr,ymax = visregUpr),alpha = .5) +
 geom_line()

The plotted figure after smoothing with RCS


编辑:您可以按天收集数据,但可以在日期中添加抖动,以使它们在一天中散布开来。

actual_df <- original_df %>%
 filter(datez %in% lubridate::date("2020-09-10")) %>%
 sample_n(size = 15) %>%
 group_by(datez) %>%
 ungroup %>%
 bind_rows(original_df %>% filter(!datez %in% lubridate::date("2020-09-10"))) %>%
 mutate(num_date = as.numeric(datez))  %>%
## Here we add random noise (uniform -.5 to .5) to each numeric date.
 mutate(jittered_date = num_date + runif(n(),-.5,.5))

## You can lower this number to increase smoothing.
number_of_knots = 15

knots <- paste0("c(",paste0(attr(rms::rcs(actual_df$jittered_date,")")

formula <- as.formula(paste("score ~ rms::rcs(jittered_date,")"))

glm_rcs <- glm(data = actual_df,family = "gaussian") %>% visreg::visreg(plot = F) %>% .$fit %>%
 mutate(date_date = chron::as.chron(jittered_date) %>% as.POSIXct())

ggplot(data = glm_rcs,y = visregFit)) +
 geom_ribbon(aes(ymin = visregLwr,alpha = .5) +
 geom_line()

figure after adding jitter to dates


编辑2:

如果您有日期时间矢量而不是简单的一天,则不必 抖动点。 在原始代码中,使用lubridate::date()创建伪造数据,它使用posix日期时间矢量并将条带化为简单日期!您可以通过以下方式避免这种情况:

original_df <- tibble(datez = seq(as.POSIXct('2020-09-01',tz = "UTC"),as.POSIXct('2020-09-15',by="1 mins") %>%
 sample(15000,replace = T)) %>%
 mutate(datez_day = lubridate::date(datez)) %>%
 add_column(score = runif(15000))

actual_df <- original_df %>%
 filter(datez_day %in% lubridate::date("2020-09-10")) %>%
 sample_n(size = 15) %>%
 bind_rows(original_df %>% filter(!datez_day %in% lubridate::date("2020-09-10"))) %>%
 mutate(num_date = as.numeric(datez))

现在您有datez_day(是日期值),datez(是日期时间)和num_date(是日期时间的数字表示)。

您可以直接在num_date上建模,而无需添加任何抖动。

number_of_knots = 20

knots <- paste0("c(",")")

formula <- as.formula(paste("score ~ rms::rcs(num_date,family = "gaussian") %>% 
        visreg::visreg(plot = F) %>% 
        .$fit %>% 
        as_tibble() %>%
   ## Translate the num_date back into a datetime object so it is correct in the figures!
        mutate(date_date = as.POSIXct.numeric(round(num_date),origin = "1970/01/01"))

ggplot(data = glm_rcs,alpha = .5) +
 geom_line()

版权声明:本文内容由互联网用户自发贡献,该文观点与技术仅代表作者本人。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如发现本站有涉嫌侵权/违法违规的内容, 请发送邮件至 dio@foxmail.com 举报,一经查实,本站将立刻删除。

相关推荐


Selenium Web驱动程序和Java。元素在(x,y)点处不可单击。其他元素将获得点击?
Python-如何使用点“。” 访问字典成员?
Java 字符串是不可变的。到底是什么意思?
Java中的“ final”关键字如何工作?(我仍然可以修改对象。)
“loop:”在Java代码中。这是什么,为什么要编译?
java.lang.ClassNotFoundException:sun.jdbc.odbc.JdbcOdbcDriver发生异常。为什么?
这是用Java进行XML解析的最佳库。
Java的PriorityQueue的内置迭代器不会以任何特定顺序遍历数据结构。为什么?
如何在Java中聆听按键时移动图像。
Java“Program to an interface”。这是什么意思?