Skip to content

climdata/globalMonthlyTI

BranchesTags

Folders and files

NameName
Last commit message
Last commit date

Latest commit

 

History

2 Commits
README_files/figure-html
README_files/figure-html
 
 
csv
csv
 
 
.gitignore
.gitignore
 
 
LICENSE
LICENSE
 
 
README.Rmd
README.Rmd
 
 
README.html
README.html
 
 
README.md
README.md
 
 

Repository files navigation

title author date output
GISS Monthly Global Temperature Index
KMicha71
17 8 2021
html_document pdf_document
keep_md
true
default

Work with yearly and moving/rolling 30y normalization

require("ggplot2")
## Loading required package: ggplot2

## Registered S3 methods overwritten by 'tibble':
##   method     from  
##   format.tbl pillar
##   print.tbl  pillar

color.temperature <- c("#0000FF", "#00CCCC", "#FFFFFF", "#EEAA33", "#FF5555")
#install.packages("data.table")
library(data.table)
#install.packages("zoo")
library(zoo)
## 
## Attaching package: 'zoo'

## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric

library(rollRegres)


tmDwd <- read.csv("https://raw.githubusercontent.com/climdata/globalMonthlyTemperature/master/csv/monthly_temperature_global.csv", sep=",")

#1949, Nov
#py3[1199,5] = 88.0 

rollingAVG <- function(data,mon,wid,col) {
  py4 <-  subset(data, (data$year>1752 & data$month == mon)) 
  avg <- rollapply(py4[col], width=wid, by=1, FUN=mean) 
  std <- rollapply(py4[col], width=wid, by=1, FUN=sd)
  py4 <-tail(py4, n=1-wid)
  py4$ind <- (py4[col] - avg)/std
  py4$ind <- signif(py4$ind, digits=6)
  return(py4)   
}


tia30 <- rollingAVG(tmDwd,1,30,'land')
for (mont in c(2,3,4,5,6,7,8,9,10,11,12)) {
  tmp <-rollingAVG(tmDwd,mont,30,'land')
  tia30 <- rbind(tia30, tmp)
}
#tia30 <- tia30[order(tia30$ts),]
tia30 <- do.call(data.frame, tia30)
names(tia30)[names(tia30) == "land.1"] <- "land.ind"

tia30o <- rollingAVG(tmDwd,1,30,'ocean')
for (mont in c(2,3,4,5,6,7,8,9,10,11,12)) {
  tmp <-rollingAVG(tmDwd,mont,30,'ocean')
  tia30o <- rbind(tia30o, tmp)
}
#tia30 <- tia30[order(tia30$ts),]
tia30o <- do.call(data.frame, tia30o)
names(tia30o)[names(tia30o) == "ocean.1"] <- "ocean.ind"

norm <- tia30
mp <- ggplot(norm, aes(year, month))
mp + geom_raster(aes(fill=land.ind))+
  scale_y_continuous(breaks=c(1,6,12))+
  theme(panel.background = element_rect(fill = '#EEEEEE', colour = 'white'), legend.position="right", text=element_text(size=14))+
  scale_fill_gradientn(colours=color.temperature)

#hist(norm$temperature)

names(tia30)[names(tia30) == "land.ind"] <- "tia30.land"
names(tia30o)[names(tia30o) == "ocean.ind"] <- "tia30.ocean"
tia30$temperature = NULL
#tia30 <- tia30[order(tia30$ts),]

Work with yearly and moving/rolling 10y normalization

tia10 <- rollingAVG(tmDwd,1,10,'land')
for (mont in c(2,3,4,5,6,7,8,9,10,11,12)) {
  tmp <-rollingAVG(tmDwd,mont,10,'land')
  tia10 <- rbind(tia10, tmp)
}
tia10 <- do.call(data.frame, tia10)
names(tia10)[names(tia10) == "land.1"] <- "land.ind"

tia10o <- rollingAVG(tmDwd,1,10,'ocean')
for (mont in c(2,3,4,5,6,7,8,9,10,11,12)) {
  tmp <-rollingAVG(tmDwd,mont,10,'ocean')
  tia10o <- rbind(tia10o, tmp)
}
tia10o <- do.call(data.frame, tia10o)
names(tia10o)[names(tia10o) == "ocean.1"] <- "ocean.ind"



norm <- tia10
mp <- ggplot(norm, aes(year, month))
mp + geom_raster(aes(fill=land.ind))+
  scale_y_continuous(breaks=c(1,6,12))+
  theme(panel.background = element_rect(fill = '#EEEEEE', colour = 'white'), legend.position="right", text=element_text(size=14))+
  scale_fill_gradientn(colours=color.temperature)

#hist(norm$temperature)

names(tia10)[names(tia10) == "land.ind"] <- "tia10.land"
names(tia10o)[names(tia10o) == "ocean.ind"] <- "tia10.ocean"
tia10$temperature = NULL
#tia10 <- tia10[order(tia10$ts),]

Work with yearly and moving/rolling 10y linearization

#install.packages("data.table")

rollingGLM <- function(data, mon, wid) {
  py4 <-  subset(data, (data$year>1752 & data$month == mon))
  reg <- roll_regres(temperature ~ year, py4, width = wid, do_compute=c('sigmas', '1_step_forecasts')) 
  #reg <- roll_regres(temperature ~ year, py4, width = wid, do_compute=c('sigmas'))
  lapply(reg, tail)
  py4$ind <- (py4$temperature - reg$one_step_forecasts)/reg$sigmas  
  py4$ind <- signif(py4$ind, digits=6)
  #py4$ind2 <- (py4$temperature - py4$year*reg$coefs[,2]+reg$coefs[,1])/reg$sigmas
  
  py4 <-tail(py4, n=-wid)
  return(py4) 
}  

tmDwd$temperature = tmDwd$land
til10 <- rollingGLM(tmDwd,1,10)
for (mont in c(2,3,4,5,6,7,8,9,10,11,12)) {
  tmp <-rollingGLM(tmDwd,mont,10)
  til10 <- rbind(til10, tmp)
}
#py5 <-tail(py5, n=-10)
til10 <- do.call(data.frame, til10)

tmDwd$temperature = tmDwd$ocean
til10o <- rollingGLM(tmDwd,1,10)
for (mont in c(2,3,4,5,6,7,8,9,10,11,12)) {
  tmp <-rollingGLM(tmDwd,mont,10)
  til10o <- rbind(til10o, tmp)
}
#py5 <-tail(py5, n=-10)
til10o <- do.call(data.frame, til10o)

norm <- til10
mp <- ggplot(norm, aes(year, month))
mp + geom_raster(aes(fill=ind))+
  scale_y_continuous(breaks=c(1,6,12))+
  theme(panel.background = element_rect(fill = '#EEEEEE', colour = 'white'), legend.position="right", text=element_text(size=14))+
  scale_fill_gradientn(colours=color.temperature)

#hist(norm$temperature)
#hist(norm$ind)

names(til10)[names(til10) == "ind"] <- "til10.land"
names(til10o)[names(til10o) == "ind"] <- "til10.ocean"
til10$temperature = NULL
#til10 <- til10[order(til10$ts),]

Work with yearly and moving/rolling 30y linearization

tmDwd$temperature = tmDwd$land
til30 <- rollingGLM(tmDwd,1,30)
for (mont in c(2,3,4,5,6,7,8,9,10,11,12)) {
  tmp <-rollingGLM(tmDwd,mont,30)
  til30 <- rbind(til30, tmp)
}

tmDwd$temperature = tmDwd$ocean
til30o <- rollingGLM(tmDwd,1,30)
for (mont in c(2,3,4,5,6,7,8,9,10,11,12)) {
  tmp <-rollingGLM(tmDwd,mont,30)
  til30o <- rbind(til30o, tmp)
}

#py5 <-tail(py5, n=-10)
norm <- til30

mp <- ggplot(norm, aes(year, month))
mp + geom_raster(aes(fill=ind))+
  scale_y_continuous(breaks=c(1,6,12))+
  theme(panel.background = element_rect(fill = '#EEEEEE', colour = 'white'), legend.position="right", text=element_text(size=14))+
  scale_fill_gradientn(colours=color.temperature)

#hist(norm$temperature)
#hist(norm$ind)

names(til30)[names(til30) == "ind"] <- "til30.land"
names(til30o)[names(til30o) == "ind"] <- "til30.ocean"
til30$temperature = NULL
#til30 <- til30[order(til30$ts),]

Work with yearly and fixed normalization

tmDwd$temperature = tmDwd$land
tis99 <-  subset(tmDwd, (tmDwd$year>1752 & tmDwd$month == 1)) 

avg <- mean(tis99$temperature, na.rm = TRUE)
std <- sd(tis99$temperature, na.rm = TRUE)
tis99$ind <- (tis99$temperature-avg)/std

for (mont in c(2,3,4,5,6,7,8,9,10,11,12)) {
  tmp <-subset(tmDwd, (tmDwd$year>1752 & tmDwd$month == mont)) 
  avg <- mean(tmp$temperature, na.rm = TRUE)
  std <- sd(tmp$temperature, na.rm = TRUE)
  tmp$ind <- (tmp$temperature-avg)/std
  tis99 <- rbind(tis99, tmp)
}
tis99$ind <- signif(tis99$ind, digits=6)


tmDwd$temperature = tmDwd$ocean
tis99o <-  subset(tmDwd, (tmDwd$year>1752 & tmDwd$month == 1)) 

avg <- mean(tis99o$temperature, na.rm = TRUE)
std <- sd(tis99o$temperature, na.rm = TRUE)
tis99o$ind <- (tis99o$temperature-avg)/std

for (mont in c(2,3,4,5,6,7,8,9,10,11,12)) {
  tmp <-subset(tmDwd, (tmDwd$year>1752 & tmDwd$month == mont)) 
  avg <- mean(tmp$temperature, na.rm = TRUE)
  std <- sd(tmp$temperature, na.rm = TRUE)
  tmp$ind <- (tmp$temperature-avg)/std
  tis99o <- rbind(tis99o, tmp)
}
tis99o$ind <- signif(tis99o$ind, digits=6)

norm <- tis99

mp <- ggplot(norm, aes(year, month))
mp + geom_raster(aes(fill=ind))+
  scale_y_continuous(breaks=c(1,6,12))+
  theme(panel.background = element_rect(fill = '#EEEEEE', colour = 'white'), legend.position="right", text=element_text(size=14))+
  scale_fill_gradientn(colours=color.temperature)

#hist(norm$temperature)

names(tis99)[names(tis99) == "ind"] <- "tis99.land"
names(tis99o)[names(tis99o) == "ind"] <- "tis99.ocean"
tis99$temperature = NULL
#tis99 <- tis99[order(tis99$ts),]

Combine indices and store

tis99 <- tis99[, c("year","month","tis99.land")]
tia30 <- tia30[, c("year","month","tia30.land")]
tia10 <- tia10[, c("year","month","tia10.land")]
til30 <- til30[, c("year","month","til30.land")]
til10 <- til10[, c("year","month","til10.land")]

tis99o <- tis99o[, c("year","month","tis99.ocean")]
tia30o <- tia30o[, c("year","month","tia30.ocean")]
tia10o <- tia10o[, c("year","month","tia10.ocean")]
til30o <- til30o[, c("year","month","til30.ocean")]
til10o <- til10o[, c("year","month","til10.ocean")]

tiAll = tis99
tiAll <- merge(tiAll,tia30, by=c("year","month"))
tiAll <- merge(tiAll,tia10, by=c("year","month"))
tiAll <- merge(tiAll,til30, by=c("year","month"))
tiAll <- merge(tiAll,til10, by=c("year","month"))

tiAll <- merge(tiAll,tis99o, by=c("year","month"))
tiAll <- merge(tiAll,tia30o, by=c("year","month"))
tiAll <- merge(tiAll,tia10o, by=c("year","month"))
tiAll <- merge(tiAll,til30o, by=c("year","month"))
tiAll <- merge(tiAll,til10o, by=c("year","month"))

tiAll$ts <- signif(tiAll$year + (tiAll$month-0.5)/12, digits=6)
tiAll$time <- paste(tiAll$year,tiAll$month, '15 00:00:00', sep='-')
tiAll <- tiAll[order(tiAll$ts),]


write.table(tiAll, file = "csv/ti_global.csv", append = FALSE, quote = TRUE, sep = ",",
            eol = "\n", na = "NA", dec = ".", row.names = FALSE,
            col.names = TRUE, qmethod = "escape", fileEncoding = "UTF-8")

Calculate multi-monthly indices and store

tiMulti = tis99
names(tiMulti)[names(tiMulti) == "tis99.land"] <- "sti1"

tiMulti$ts <- signif(tiMulti$year + (tiMulti$month-0.5)/12, digits=6)
tiMulti$time <- paste(tiMulti$year,tiMulti$month, '15 00:00:00', sep='-')
tiMulti <- tiMulti[order(tiMulti$ts),]

prev <- tiMulti$sti1
for (m in c(2,3,4,5,6,7,8,9,10,11,12)) {
  column <- paste("sti", m, sep="")
  sti <- rollapply(tiMulti$sti1, width=m, by=1, FUN=sum)
  tiMulti$sti <- prev
  tiMulti$sti[m:length(tiMulti$sti)] <- sti
  tiMulti$sti <- tiMulti$sti*m^(-1/sqrt(3))
  prev <- tiMulti$sti  
  names(tiMulti)[names(tiMulti) == 'sti'] <- column
}
tiMulti <- tiMulti[order(tiMulti$ts),]

write.table(tiMulti, file = "csv/sti_de.csv", append = FALSE, quote = TRUE, sep = ",",
            eol = "\n", na = "NA", dec = ".", row.names = FALSE,
            col.names = TRUE, qmethod = "escape", fileEncoding = "UTF-8")

Plot monthly

require("ggplot2")
temp <- read.csv("./csv/ti_global.csv", sep=",")

mp <- ggplot() +
      geom_line(aes(y=temp$tis99.ocean, x=temp$ts), color="blue") +
      geom_line(aes(y=temp$tis99.land, x=temp$ts), color="brown") +  
      xlab("Year") + ylab("tis99 ")
mp

mp <- ggplot() +
      geom_line(aes(y=temp$tia30.ocean, x=temp$ts), color="blue") +
      geom_line(aes(y=temp$tia30.land, x=temp$ts), color="brown") +  
      xlab("Year") + ylab("tia30 ")
mp

mp <- ggplot() +
      geom_line(aes(y=temp$tia10.ocean, x=temp$ts), color="blue") +
      geom_line(aes(y=temp$tia10.land, x=temp$ts), color="brown") +  
      xlab("Year") + ylab("tia10 ")
mp

mp <- ggplot() +
      geom_line(aes(y=temp$til30.ocean, x=temp$ts), color="blue") +
      geom_line(aes(y=temp$til30.land, x=temp$ts), color="brown") +  
      xlab("Year") + ylab("til30 ")
mp

mp <- ggplot() +
      geom_line(aes(y=temp$til10.ocean, x=temp$ts), color="blue") +
      geom_line(aes(y=temp$til10.land, x=temp$ts), color="brown") +  
      xlab("Year") + ylab("til10 ")
mp

xxx

require("ggplot2")
temp <- read.csv("./csv/ti_global.csv", sep=",")

mp <- ggplot() +
      geom_line(aes(y=temp$til30.land, x=temp$ts), color="green") +  
      geom_line(aes(y=temp$tia10.land, x=temp$ts), color="cyan") +  
      geom_line(aes(y=temp$til10.land, x=temp$ts), color="magenta") +  
      geom_line(aes(y=temp$tia30.land, x=temp$ts), color="red") +  
      geom_line(aes(y=temp$tis99.land, x=temp$ts), color="blue") +
      xlab("Year") + ylab("land ")
mp

mp <- ggplot() +
      geom_line(aes(y=temp$til30.ocean, x=temp$ts), color="green") +  
      geom_line(aes(y=temp$tia10.ocean, x=temp$ts), color="cyan") +  
      geom_line(aes(y=temp$til10.ocean, x=temp$ts), color="magenta") +  
      geom_line(aes(y=temp$tia30.ocean, x=temp$ts), color="red") +  
      geom_line(aes(y=temp$tis99.ocean, x=temp$ts), color="blue") +
      xlab("Year") + ylab("ocean ")
mp

About

No description, website, or topics provided.

Resources

Readme

License

Apache-2.0 license
Activity
Custom properties

Stars

0 stars

Watchers

1 watching

Forks

0 forks
Report repository

Releases

No releases published

Packages

No packages published