Climate_change_all_Densities: CB14 and CB15 move

setwd("D:/Compositional regression/July2022")
load("HclimateB2.RDATA")
require(tidyverse)
require(graphics)

#=============Plot the results
yrname <- "2016"
provincename <- unique(names(HclimateB2[[yrname]]))[50]#[1:10]

spline0B2_n <- HclimateB2[[yrname]][[provincename]]

n <- 2
h12_n <- c(1, -1 )

Year 2014

#op <- par(mfcol=c(1,1))

yrname <- "2014"
  for (provincename in unique(names(HclimateB2[[yrname]])))
       {
         spline0B2_n <- HclimateB2[[yrname]][[provincename]]
         matplot(spline0B2_n$Partition, 
        spline0B2_n[, 3:(n+2)], 
        type = "l", lty = 1, #las = 1, 
        ylab = "Density", xlab = "Temperature",
        lwd = 1.5, cex.lab = 1, 
        cex.axis = 0.5, 
        col = rainbow(7)[1:n],#)#,
   main = c(paste("Province", provincename , "and year", yrname ))) #and 16th
lines(spline0B2_n[, 1], 
      spline0B2_n[, 2], lwd = 1.5)
#abline(v = CFR.tmax$knots, col = "gray", lty = 2)
legend("topleft", legend = c(paste0( "h = ", h12_n), "Original"), 
       col =  c(rainbow(7)[1:n], "black"),
       ncol = 2, cex = 0.5, lwd = 1.5)
  }

#par(op)

Year 2015

#op <- par(mfcol=c(1,1))

yrname <- "2015"
  for (provincename in unique(names(HclimateB2[[yrname]])))
       {
         spline0B2_n <- HclimateB2[[yrname]][[provincename]]
         matplot(spline0B2_n$Partition, 
        spline0B2_n[, 3:(n+2)], 
        type = "l", lty = 1, #las = 1, 
        ylab = "Density", xlab = "Temperature",
        lwd = 1.5, cex.lab = 1, 
        cex.axis = 0.5, 
        col = rainbow(7)[1:n],#)#,
   main = c(paste("Province", provincename , "and year", yrname ))) #and 16th
lines(spline0B2_n[, 1], 
      spline0B2_n[, 2], lwd = 1.5)
#abline(v = CFR.tmax$knots, col = "gray", lty = 2)
legend("topleft", legend = c(paste0( "h = ", h12_n), "Original"), 
       col =  c(rainbow(7)[1:n], "black"),
       ncol = 2, cex = 0.5, lwd = 1.5)
  }

#par(op)

Year 2016

#op <- par(mfcol=c(1,1))

yrname <- "2016"
  for (provincename in unique(names(HclimateB2[[yrname]])))
       {
         spline0B2_n <- HclimateB2[[yrname]][[provincename]]
         matplot(spline0B2_n$Partition, 
        spline0B2_n[, 3:(n+2)], 
        type = "l", lty = 1, #las = 1, 
        ylab = "Density", xlab = "Temperature",
        lwd = 1.5, cex.lab = 1, 
        cex.axis = 0.5, 
        col = rainbow(7)[1:n],#)#,
   main = c(paste("Province", provincename , "and year", yrname ))) #and 16th
lines(spline0B2_n[, 1], 
      spline0B2_n[, 2], lwd = 1.5)
#abline(v = CFR.tmax$knots, col = "gray", lty = 2)
legend("topleft", legend = c(paste0( "h = ", h12_n), "Original"), 
       col =  c(rainbow(7)[1:n], "black"),
       ncol = 2, cex = 0.5, lwd = 1.5)
  }

#par(op)