ProgramBox

Как рассчитать отклонение температуры от нормы, используя среднее значение за 30 лет — R 4.0.0

Post author:admin
Запись опубликована:13 марта, 2023
Post category:Вопросы по программированию

#r #rnoaa

Вопрос:

Я добавил весь свой код ниже, так что я проверил, работает ли он:

Я точно знаю, что могу загрузить любую климатическую станцию с веб-сайта NCDC для определенного диапазона времени, который я хочу. Кстати, если вы можете взглянуть на мою команду ‘bind_rows ()` и сделать ее менее запутанной, я не смог найти лучшего способа сделать это.
Я знаю TAVG , что вычисляется и работает
Ежемесячные сводки, в результате чего data.set mso_sum работает отлично

Так что же не работает:

Нахождение моих отклонений от 30-летних норм

Как я хочу, чтобы это работало:

отфильтруйте годы 1981: 2010
группируйте так, чтобы каждое 1 января, 2 января, 3 января и т. Д., 1 февраля, 2 февраля и т. Д. Можно суммировать по дням
суммируйте среднее значение TAVG (среднее значение температуры, полученное из maxT и MinT)
затем возьмите весь набор данных и вычтите ежедневный TAVG из CliAvgT

Это код, который я пробовал:

 mso_light %>%
  group_by(month, day) %>%
  summarise(CliAvgT = mean(TAVG[1981:2010], na.rm = T)) %>%
  mutate(Avg_DepT = CliAvgT - TAVG) %>%
  ungroup()

Я также попробовал этот альтернативный код:

 mso_light %>%
  filter(year >= "1981", year <= "2010") %>%
  group_by(month, day) %>%
  summarise(CliAvgT = mean(TAVG, na.rm = T)) %>%
  mutate(Avg_DepT = CliAvgT - TAVG) %>%
  ungroup()

и получил это сообщение об ошибке:

 > mso_light %>%
Warning messages:
1: Unknown or uninitialised column: `value`. 
2: Unknown or uninitialised column: `value`. 
    filter(year >= "1981", year <= "2010") %>%
    group_by(month, day) %>%
    summarise(CliAvgT = mean(TAVG, na.rm = T)) %>%
    mutate(Avg_DepT = CliAvgT - TAVG) %>%
    ungroup()
`summarise()` regrouping output by 'month' (override with `.groups` argument)
Error: Problem with `mutate()` input `Avg_DepT`.
x object 'TAVG' not found
i Input `Avg_DepT` is `CliAvgT - TAVG`.
i The error occured in group 1: month = 1.
Run `rlang::last_error()` to see where the error occurred.

Наконец, вот весь мой код:

 library(rnoaa)
library(tidyverse)
library(data.table)
library("openair")
library("chron")
library('lubridate')

## grab first half of year

getNoaaP1 <- function(yr, type = c('tmax','tmin','PRCP', 'SNOW', 'SNWD')) 
  ncdc(datasetid = 'GHCND', 
       stationid = 'GHCND:USW00024153',
       datatypeid = type, 
       startdate = paste0(yr, '-01-01'), 
       enddate = paste0(yr, '-06-30'), limit = 1000)  

## grab second half of year

getNoaaP2 <- function(yr, type = c('tmax','tmin','PRCP', 'SNOW', 'SNWD')) 
  ncdc(datasetid = 'GHCND', 
       stationid = 'GHCND:USW00024153',
       datatypeid = type, 
       startdate = paste0(yr, '-07-01'), 
       enddate = paste0(yr, '-12-31'), limit = 1000)  

res1 <- setNames(lapply(1948:2020, getNoaaP1), paste0("Year", 1948:2020, "P1"))
res <- setNames(lapply(1948:2020, getNoaaP2), paste0("Year", 1948:2020, "P2"))

# this would export all individual list elements to the global environment:
list2env(res, envir = .GlobalEnv) 
list2env(res1, envir = .GlobalEnv) 

# this would combine the individual lists
mso <- bind_rows(Year1948P1$data, Year1949P1$data, Year1950P1$data, Year1951P1$data,
                 Year1952P1$data, Year1953P1$data, Year1954P1$data, Year1955P1$data,
                 Year1956P1$data, Year1957P1$data, Year1958P1$data, Year1959P1$data,
                 Year1960P1$data, Year1961P1$data, Year1962P1$data, Year1963P1$data,
                 Year1964P1$data, Year1965P1$data, Year1966P1$data, Year1967P1$data,
                 Year1968P1$data, Year1969P1$data, Year1970P1$data, Year1971P1$data,
                 Year1972P1$data, Year1973P1$data, Year1974P1$data, Year1975P1$data,
                 Year1976P1$data, Year1977P1$data, Year1978P1$data, Year1979P1$data,
                 Year1980P1$data, Year1981P1$data, Year1982P1$data, Year1983P1$data,
                 Year1984P1$data, Year1985P1$data, Year1986P1$data, Year1987P1$data,
                 Year1988P1$data, Year1989P1$data, Year1990P1$data, Year1991P1$data,
                 Year1992P1$data, Year1993P1$data, Year1994P1$data, Year1995P1$data,
                 Year1996P1$data, Year1997P1$data, Year1998P1$data, Year1999P1$data,
                 Year2000P1$data, Year2001P1$data, Year2002P1$data, Year2003P1$data,
                 Year2004P1$data, Year2005P1$data, Year2006P1$data, Year2007P1$data,
                 Year2008P1$data, Year2009P1$data, Year2010P1$data, Year2011P1$data,
                 Year2012P1$data, Year2013P1$data, Year2014P1$data, Year2015P1$data,
                 Year2016P1$data, Year2017P1$data, Year2018P1$data, Year2019P1$data,
                 Year2020P1$data,
                 Year1948P2$data, Year1949P2$data, Year1950P2$data, Year1951P2$data,
                 Year1952P2$data, Year1953P2$data, Year1954P2$data, Year1955P2$data,
                 Year1956P2$data, Year1957P2$data, Year1958P2$data, Year1959P2$data,
                 Year1960P2$data, Year1961P2$data, Year1962P2$data, Year1963P2$data,
                 Year1964P2$data, Year1965P2$data, Year1966P2$data, Year1967P2$data,
                 Year1968P2$data, Year1969P2$data, Year1970P2$data, Year1971P2$data,
                 Year1972P2$data, Year1973P2$data, Year1974P2$data, Year1975P2$data,
                 Year1976P2$data, Year1977P2$data, Year1978P2$data, Year1979P2$data,
                 Year1980P2$data, Year1981P2$data, Year1982P2$data, Year1983P2$data,
                 Year1984P2$data, Year1985P2$data, Year1986P2$data, Year1987P2$data,
                 Year1988P2$data, Year1989P2$data, Year1990P2$data, Year1991P2$data,
                 Year1992P2$data, Year1993P2$data, Year1994P2$data, Year1995P2$data,
                 Year1996P2$data, Year1997P2$data, Year1998P2$data, Year1999P2$data,
                 Year2000P2$data, Year2001P2$data, Year2002P2$data, Year2003P2$data,
                 Year2004P2$data, Year2005P2$data, Year2006P2$data, Year2007P2$data,
                 Year2008P2$data, Year2009P2$data, Year2010P2$data, Year2011P2$data,
                 Year2012P2$data, Year2013P2$data, Year2014P2$data, Year2015P2$data,
                 Year2016P2$data, Year2017P2$data, Year2018P2$data, Year2019P2$data,
                 Year2020P2$data)

## build data.frame and remove 'station ID' column
mso_light <- mso[, -3]

## remove time from date group
mso_date <- mso_light[1]
mso_date <- sub("T.*", "", mso_date$date)
mso_light$date <- mso_date 

## remove flags for fl_so? and fl_t (time)
mso_light <- mso_light[1:5]

## Change 'T' = 9998 amp; 'M' = 9999
mso_light$value[mso_light$fl_m == "T"] <- 0
mso_light$value[mso_light$fl_q == "M"] <- 'na'

mso_light$value <- as.numeric(mso_light$value)

## pivot data frame

## eventually use to change column names
## v_names <- c('PRCP', 'SNOW', 'SNWD', 'TMAX', 'TMIN')

mso_light <- mso_light[1:3]

mso_light <- pivot_wider(mso_light,
                          names_from = datatype,  
                          values_from = value)
## mso_light <- select(mso_light, -c('fl_m','fl_q'))

options(stringAsFactors = FALSE)

mso_light$date <- as.Date(mso_light$date, "%Y-%m-%d")


## Turning all daily temperatures into an average

mso_light <- mso_light %>% rowwise() %>%
  mutate(TAVG = mean(c(TMAX, TMIN), na.rm = T))

## Composing daily data into monthly packages

mso_light <- mso_light %>%
  mutate(month = month(date)) %>%
  mutate(year = year(date)) %>%
  mutate(day = day(date))

mso_light <- mso_light %>%
  relocate('year', 'month', 'day') 

## mso_light <- mso_light[-4]

mso_sum <- mso_light %>%
  group_by(month, year) %>% 
  summarize(AVG_TAVG=mean(TAVG, na.rm = TRUE),
          T_PRCP=sum(PRCP, na.rm=TRUE),
          T_SNOW=sum(SNOW, na.rm=TRUE)) %>% 
  ungroup()

## make 30 year averages, using 1981-2010

mso_light %>%
  group_by(month, day) %>%
  summarise(CliAvgT = mean(TAVG[1981:2010], na.rm = T)) %>%
  mutate(Avg_DepT = CliAvgT - TAVG) %>%
  ungroup()

##mso_DeptT <- mso_light %>%
##  group_by(month, day) %>%
##  mean(mso_light$TAVG[1981:2010], na.rm = T) %>%
##  ungroup()

##mso_DeptT <- filter(mso_light, year >= "1981", year <= "2010") %>%
##  group_by(day, month) %>%
##  mutate(daily_DeptT = mean(TAVG, na.rm = T)) %>%
##  ungroup()

cli_Avg <- filter(mso_sum, year >= "1981", year <= "2010") %>%
  group_by(month) %>%
  summarize(T_dep = mean(AVG_TAVG, na.rm = T),
            Mon_Precip = mean(T_PRCP, na.rm = T),
            Mon_Snow = mean(T_SNOW, na.rm = T))

write.csv(mso_light, "mso_light.csv")
write.csv(mso_sum, "mso_sum.csv")
write.csv(cli_Avg, "cli_avg.csv")

Ответ №1:

использование data.table вашего кода может быть написано как таковое. Также проблема в том, что summarise удаляет все остальные столбцы, и вы могли бы использовать его следующим образом summarise(CliAvgT = mean(TAVG, na.rm = T))

аккуратный способ

 mso_light %>%
  group_by(month, day) %>%
  summarise(CliAvgT = mean(TAVG, na.rm = T)) %>%
  ungroup() %>% right_join(mso_light) %>%
  mutate(Avg_DepT = CliAvgT - TAVG)

использование `data.table` синтаксиса

это полная переформулировка вашего кода :

 # rbinding all the data into one data.table
# lapply(c(res,res1),`[[`, "data") <=> lapply(c(res,res1),function(x) x[["data"]])
mso <- rbindlist(lapply(c(res,res1),`[[`, "data"))

# remove the station column
# remove time from date 
mso[,c("station", "date"):=list(NULL, as.Date(sub("T. $", "", date)))]

# dcast <=> pivot_wider
# mso[,1:3] get the first three columns
# and pivot them wider 
# pivot_wider(names_from, values_from) <=> dcast(., idcol ~ names_from, value.var=values_from)
mso_light <- dcast(mso[,1:3], date ~ datatype, value.var = "value")

# calculate TAVG as we only have two values => TAVG= (TMAX TMIN)/2
# then create the year, month day columns
mso_light[, TAVG:= (TMAX TMIN)/2][, c("year", "month", "day") := list(year(date),month(date),day(date))]

# create a new data.table that contains the year, month and the average TAVG the sums of PRCP and SNOW by month and year
mso_sum <- mso_light[,list(year, month, AVG_TAVG=mean(TAVG, na.rm = TRUE),
          T_PRCP=sum(PRCP, na.rm=TRUE),
          T_SNOW=sum(SNOW, na.rm=TRUE)), by=c("year","month")]

# 
mso_light[, CliAvgT := mean(TAVG, na.rm = T), by=c("month", "day")][, Avg_DepT := CliAvgT - TAVG]

cli_Avg <- mso_sum[year>=1981 amp; year<=2020, list(T_dep = mean(AVG_TAVG, na.rm = T),
            Mon_Precip = mean(T_PRCP, na.rm = T),
            Mon_Snow = mean(T_SNOW, na.rm = T)), by=month ]

некоторые пояснения:

в принципе, использование data.table [ имеет несколько новых аргументов. Я постараюсь дать краткое введение :

df[, new.column:=old.col*2] : создает new.column и добавляет его к текущему df <=> df$new.column <- df$old.col*2
df[col<0] возвращает строки, которые проверяют col<0 , т.е. df[df$col<0]
вы можете смешать их, и это в основном изменит только те строки, в которых проверяется условие.
by аргумент: в основном группирует данные по предоставленным столбцам df[, , by=col]<=> df %>% group_by(col)
последнее предостережение — df[, c(col1,col2)] вернет вектор, содержащий значения в этих столбцах, а использование list вместо c вернет data.table, содержащий значения столбцов.

вводный урок

1. Единственная проблема, которую я вижу сразу, заключается в том, что «отклонение от нормы» взято из 1981-2010 годов. То, как вы это написали, охватило бы весь набор данных. Как мне отфильтровать его по этим конкретным годам 1981-2010? Мне придется идти построчно и разбивать то, что вы сделали, но это выглядит намного проще для чтения.

2. хорошо, прежде всего [[ , это просто базовая R функция подмножества, то есть: ` [[ (list.or.vec, element.index) <=> list.or.vec[[element.index]]`

3. @Metgeneer Я добавил несколько объяснений

Как рассчитать отклонение температуры от нормы, используя среднее значение за 30 лет — R 4.0.0

Вопрос:

Комментарии:

Ответ №1:

аккуратный способ

использование `data.table` синтаксиса

некоторые пояснения:

Комментарии:

Вопрос:

Комментарии:

Ответ №1:

аккуратный способ

использование data.table синтаксиса

некоторые пояснения:

Комментарии:

Вам также может понравиться

Как использовать компонент splitter в фреймворке quasar для динамического вырезания страниц?

отправка почты через gmail

Создание td-контракта только для его заполнения и содержимого

использование `data.table` синтаксиса