Snødybde
Import
dat <- read_excel("P:/41201042_okologisk_tilstand_fastlandsnorge_2020_dataanaly/fjell2021/data/Klima/Snodybde/snowDepth_med.xlsx")dat$year2 <- as.numeric(substr(dat$year, 6, 10))
summary(dat$year2)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 1958 1973 1989 1989 2005 2020dat <- dat[dat$year>1960,]dat$reg <- as.factor(dat$reg)
table(dat$reg)##
## Midt-Norge Nord-Norge Sørlandet Vestlandet Østlandet
## 60 60 60 60 60Verdiene er i mm og er regnet ut som gjennomsnitet per celle i perioden des-mai. Deretter har man tatt medianen av cellene.
head(dat$value)## [1] 626.8187 660.0055 917.5330 794.4973 753.2967 597.5275Regner ut referanseverdien og 2SD av denne.
ref <- aggregate(data =
dat[dat$year2 %between% c(1961, 1990),],
value~reg,
FUN = mean)
upp <- aggregate(data =
dat[dat$year2 %between% c(1961, 1990),],
value~reg,
FUN = sd)
upp$value <- upp$value*2
ref$upp <- upp$value
rm(upp)Trender
regOrder = c(
"Nord-Norge",
"Midt-Norge",
"Østlandet",
"Sørlandet",
"Vestlandet"
)dat$ref <- ref$value[match(dat$reg, ref$reg)]
dat$diff <- dat$value-dat$ref
dat$col <- ifelse(dat$diff<0, "one", "two")brk <- dat$year2[seq(5, 60, 5)]
lab <- dat$year[seq(5, 60, 5)] gg <- ggplot(data = dat,
aes(x = year2, y = diff))+
geom_bar(stat="identity", aes( fill = col))+
geom_hline(data = ref,
aes(yintercept = -upp),
linetype=2)+
geom_hline(yintercept = 0)+
geom_smooth(data = dat,
aes(x = year2, y = diff))+
scale_fill_hue(l=70, c=60)+
theme_bw(base_size = 20)+
ylab("Snødybde (mm)\n avvik fra 1961-1990")+
xlab("")+
theme(
axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1))+
scale_x_continuous(breaks = brk, labels = lab)+
guides(fill="none")+
facet_wrap( .~ factor(reg, levels = regOrder),
ncol=3)png("../output/indicatorPlots/supporting figures/snødybde normalisert tiddserie med barplot.png",
units="in", width=5, height=12, res=300)
gg
dev.off()
Stiplalinja er 2sd under forrige normalperioden.
Her er en forenklet figur i ØT-stilen.
dat$yearLabel <- dat$year
dat$year <- dat$year2regOrder = c("Østlandet","Sørlandet","Vestlandet","Midt-Norge","Nord-Norge")
dat <- dat[order(match(dat$reg,regOrder),dat$year),]
minyear <- 1958
maxyear <- 2021
upperYlimit <- 1000
lowYlimit <- -1000
yStep <- 200
move <- 0.2
legendPosition <- "top"
legendInset = 0
horizontal = TRUE
legendTextSize = 1.25
colours = c("#2DCCD3", "#004F71", "#7A9A01", "#93328E", "#FFB25B")
# Create loop factors
uniq1 <- unique(unlist(dat$year))
uniq2 <- unique(unlist(dat$reg))
### PLOT first Norway
# Subset for region 'Norge'
Norge <- subset(dat, reg=="Østlandet")
png("../output/indicatorPlots/uskalert/snødybde.png",
units="in", width=12, height=7, res=300)
par(mar=c(4.5,6.5,2,2))
# Plot for region = 'Norge'
plot(
Norge$diff~Norge$year,
ylab="Snødybde (mm)\n avvik fra normalperiodn (1961-1990)",
xlab="",
main="",
xlim=c(minyear, maxyear),
ylim=c(lowYlimit, upperYlimit),
cex.main=1,
cex.lab=1.5,
cex.axis=1.5,
type="n",
frame.plot=FALSE,
axes=FALSE
)
# Axis 1 options
axis(side=1, at=seq(1960, 2020, by=10), cex.axis=1.5)
# Axis 2 options
axis(side=2, at=seq(lowYlimit, upperYlimit, yStep),
labels=seq(lowYlimit, upperYlimit, yStep),
cex.axis=1.5)
# Add lines
lines(Norge$year+(move*(-2.5)), Norge$diff, col=colours[5], lwd=2, lty=1)
# Save temp points for later addition to plot
temppoints <- data.frame(year = Norge$year, med = Norge$diff)
# Empty temporary points data frame
temppoints3 <- data.frame()
### Then plot loop per region
for(n in 1:(length(uniq2)-1)){
# Subset for region i
quants <- subset(dat, reg==uniq2[n])
# Add lines
lines(quants$year+move*(n-2.5), quants$diff, col=colours[n], lwd=2, lty=1)
# Save temp points for later addition to plot
temppoints2 <- data.frame(year = quants$year, med = quants$diff, reg = uniq2[n])
temppoints3 <- rbind(temppoints3, temppoints2)
}
## Add points for regions
#for(n in 1:(length(uniq2)-1)){
# temppoints4 <- temppoints3[temppoints3$reg==uniq2[n],]
# points(temppoints4$year+move*(n-2.5),temppoints4$diff, pch=21, bg=colours[n], cex=1.5)
#}
#
## Add points for Norge
#points(temppoints$year+(move*(-2.5)),temppoints$diff, pch=21, bg=colours[6], cex=1.5)
# Add legend to plot
legend(legendPosition, legendPositionY, legend = regOrder, col = c(colours[5], colours[1:4]),
#bg = c(colours),
pch=16, lty=2,
lwd=1.5, bty="n", inset=legendInset, title="", horiz = horizontal,
cex=legendTextSize)
# add reference line
abline(h=0, col="black", lwd=2, lty=2)
dev.off()
Bootstrapping (60% uten tilbakelegging) av gjennomsnittlig vinterregn siste 5 år.
new <- dat[dat$year2 %between% c(2016, 2020),]
round(tapply(new$diff, new$reg, FUN = mean), 1)## Midt-Norge Nord-Norge Sørlandet Vestlandet Østlandet
## -24.2 158.4 -108.8 -20.5 -27.9Tabellen over viser gjennomsnittlig reduksjon i snødybde siste 5 år.
sno <- data.frame(
reg = rep(levels(dat$reg), each = 10000),
year = 2020,
val = NA)
for(n in levels(dat$reg)){
temp <- new[new$reg==n,]
for(i in 1:10000){
sno$val[i+10000*(which(levels(dat$reg)==n)-1)] <-
mean(sample(temp$value, 3, replace=F))
}
}Skalering
Jeg setter den min-verdi som null dager. Dette betyr at indikatoren får verdi null når det ikke lenger er noe vinter i fjellet.
sno$ref <- ref$value[match(sno$reg, ref$reg)]
sno$err <- ref$upp[match(sno$reg, ref$reg)]
sno$extr <- sno$ref-sno$err
# trunkerer
sno$valS <- ifelse( sno$val>sno$ref,
sno$ref,
sno$val)
any(sno$valS<0)## [1] FALSElow <- sno[sno$valS<sno$extr,]
upp <- sno[sno$valS>=sno$extr,]
upp$val_scaled <- 1-((upp$ref-upp$valS)/(upp$ref-upp$extr))
upp$val_scaled2 <- 0.6+(upp$val_scaled*0.4)
low$val_scaled <- 1-((low$extr-low$valS)/(low$extr))
low$val_scaled2 <- low$val_scaled*0.6
sno_s <- rbind(upp, low)
sno_s <- select(sno_s, reg, year, val = val_scaled2)ggplot(data=sno_s, aes(x = val))+
geom_histogram()+
facet_wrap(.~reg)## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
Arealvekting
wgt <- readRDS("../data/fjellareal.rds")
wgt$Fjellareal2 <- wgt$Fjellareal/max(wgt$Fjellareal)
wgt$reg <- c("N", "C", "E", "W", "S")norge <- data.frame(
reg = rep("Norge", 10000),
year = 2020,
val = NA
)
temp <- sno_s
temp2 <- c(
sample(temp$val[temp$reg == "Nord-Norge"], wgt$Fjellareal2[wgt$reg == "N"]*10000, replace =T),
sample(temp$val[temp$reg == "Østlandet"], wgt$Fjellareal2[wgt$reg == "E"]*10000, replace =T),
sample(temp$val[temp$reg == "Vestlandet"], wgt$Fjellareal2[wgt$reg == "W"]*10000, replace =T),
sample(temp$val[temp$reg == "Sørlandet"], wgt$Fjellareal2[wgt$reg == "S"]*10000, replace =T),
sample(temp$val[temp$reg == "Midt-Norge"], wgt$Fjellareal2[wgt$reg == "C"]*10000, replace =T)
)
temp3 <- sample(temp2, 10000, replace = F)
norge$val <- temp3sno2 <- rbind(sno_s, norge)Plotting
myPlot <- ggplot(data = sno2,
aes(x = factor(reg, levels = regOrder),
y = val))+
geom_boxplot(fill = "grey", lwd=1.2)+
ylab("Gjennomsnittlig snødybde\nskalert mot referanseverdi")+
xlab("")+
scale_y_continuous(breaks = seq(0, 1, 0.2), limits=c(0,1))+
geom_hline(yintercept = 0.6, size=1.2, linetype=2)+
theme_bw(base_size = 20)+
theme(
axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1)
)
png("../output/indicatorPlots/skalert/snødybde_boxplot.png",
units="in", width=5, height=7, res=300)
myPlot
dev.off()
Tabell
sno2$reg <- revalue(sno2$reg,
c("Nord-Norge"="N",
"Midt-Norge"="C",
"Østlandet" ="E",
"Vestlandet"="W",
"Sørlandet"="S"))
unique(sno2$reg)## [1] "C" "N" "S" "W" "E" "Norge"sno2$X <- "snodybde"Export csv
write.csv(sno2, "../output/indicator_values/snodybde.csv", row.names = F)