美文网首页机器学习与计算机视觉R语言与统计分析数据-R语言-图表-决策-Linux-Python
R - 核支持向量机预测模型构建详解:以短波辐射通量shortw

R - 核支持向量机预测模型构建详解:以短波辐射通量shortw

作者: TroyShen | 来源:发表于2019-12-30 17:25 被阅读0次

    目录

    • 0.问题导入
    • 1.示例数据
    • 2.导入示例数据
    • 3.将示例数据按照70%/30%拆分为训练组与验证组
    • 4.构建核支持向量机预测模型(SVR)
    • 5.模型训练结果与实测数据相关性及RMSE分析(图1)
    • 6.模型训练结果与实测数据时间路径分析(图2)
    • 7.SVR模型预测结果与验证集实测相关性及RMSE分析(图3)
    • 8.SVR模型预测结果与验证集实测数据时间路径分析(图4)
    • 9.模型模拟残差分析(图5)
    • 10.总结
    • 11.本篇所使用的软件包(没有的需要通过install.packages('包名')进行安装)
    • 12.致谢

    0. 问题导入

    今天我们通过建立短波辐射通量与温度及蒸散发的关系来详细阐述核支持向量机预测模型(SVR)在面向地理及时间序列数据建模方面的应用。

    1.示例数据

    本篇数据与R - 广义加性模型(GAM)构建详解:以归一化植被指数NDVI为例这篇文章数据一致,时间长度为2002-01 ~ 2015-09,时间分辨率为月。
    点我下载示例数据
    示例数据预览:

    head(data)
              ndvi       soil1      soil2     soil3      soil4        rain
    1  0.237984169 -0.01210715 0.03579731 0.1269299 0.07318894 -0.01543584
    2  0.370455335  0.38147139 0.31089661 0.2241396 0.10204067  0.20701857
    3  0.331657733  0.41044975 0.48385978 0.4471074 0.25112199  0.62105802
    4  0.216662956  0.32583872 0.41198999 0.4231082 0.42613716  0.37216417
    5  0.054132382  0.24177292 0.20540345 0.2979310 0.43549429  0.06553887
    6 -0.005636952  0.41268755 0.29207486 0.2508858 0.37087816  0.25502620
          longwave    shortwave    root_sm        evap      temper         gpp
    1  0.059414987  0.215758745 0.06890271 -0.07747205 0.009909431 -0.04072053
    2  0.009142641  0.244385277 0.31129426  0.23793998 0.172678808  0.18329118
    3 -0.097150000  0.353491078 0.48706357  0.59985033 0.314583437  0.24478460
    4 -0.031527285  0.355970841 0.42948289  0.51469995 0.348457057  0.47172457
    5 -0.291598633  0.297255464 0.27643746  0.38420614 0.326270291  0.37802032
    6  0.010729154 -0.009709685 0.32028271  0.31684306 0.119881673  0.15986512
    

    2. 导入示例数据

    setwd('/Users/jerseyshen/Documents/JianShu_Project/20191230')
      data = read.csv("data.csv",
                      header = T)
      data = data[,-1]
      date = seq(as.Date('2002-01-01'),
                 as.Date('2015-12-01'),
                 '1 month')[1:165]
    

    3. 将示例数据按照70%/30%拆分为训练组与验证组

     train = data[1:115,]
      test = data[116:165,]
      date_train = date[1:115]
      date_test = date[116:165]
    

    4. 构建核支持向量机预测模型(SVR)

    本例采用的SVR模型中,回归模型为:eps-regression ; kennel为:RBF(Radial Base Function)。

    model_init = svm(shortwave ~ temper+evap,
                       data = train,
                       kernel = 'radial')
    

    5. 模型训练结果与实测数据相关性及RMSE分析(图1)

    train_simu = as.numeric(model_init$fitted)
      pl_point_df = data.frame(SIMU = as.numeric(train_simu),
                               TRAIN = train$shortwave)
      label_df = data.frame(
        x = c(-0.5,-0.5,-0.5),
        y = c(0.3,0.25,0.2),
        label = c(paste0("COR = ", round(cor(train_simu,train$shortwave),2)),
                  paste0("p = ", formatC(cor.test(train_simu,train$shortwave)$p.value)),
                  paste0("RMSE = ", round(rmse(train_simu,train$shortwave),3)))
      )
      
      p1_cor = ggplot()+
        geom_point(data = pl_point_df,
                   aes(x = SIMU, y= TRAIN),
                   size = 5,
                   color = 'blue',
                   shape = 1)+
        geom_abline(intercept = 0,slope = 1,size = 1)+
        geom_text(data = label_df,aes(x = x,y = y,label = label),size = 6,color = 'black',
                  hjust = 0)+
        theme_bw()+
        theme(
          axis.text =  element_text(face = 'bold',color = 'black',size = 12),
          axis.title =  element_text(face = 'bold',color = 'black',size = 14, hjust = .5)
        )+
        xlab('SIMULATION RESULT')+
        ylab('REAL SHORTWAVE')
      
      png('plot1.png',
          height = 20,
          width = 20,
          units = 'cm',
          res = 800)
      print(p1_cor)
      dev.off()
    
    图1 模型训练结果与实测数据相关性及RMSE分析

    6. 模型训练结果与实测数据时间路径分析(图2)

    pl_line_df = data.frame(
        Date = date_train,
        SIMU = train_simu,
        REAL = train$shortwave
      )
      pl_line_df = melt(pl_line_df,'Date')
      
      p1_line = ggplot()+
        geom_line(data = pl_line_df,
                   aes(x = Date, y= value,color= variable,linetype = variable),
                   size= 1)+
        scale_color_manual(values = c('green','blue'))+
        scale_linetype_manual(values = c('solid','dashed'))+
        theme_bw()+
        theme(
          axis.text =  element_text(face = 'bold',color = 'black',size = 12),
          axis.title =  element_text(face = 'bold',color = 'black',size = 14, hjust = .5),
          legend.text = element_text(face = 'bold',color = 'black',size = 12),
          legend.title = element_blank(),
          legend.position = 'bottom',
          legend.direction = 'horizontal'
        )+
        xlab('Date')+
        ylab('Indices')
      
      png('plot2.png',
          height = 10,
          width = 20,
          units = 'cm',
          res = 800)
      print(p1_line)
      dev.off()
    
    图2 模型训练结果与实测数据时间路径分析

    7. SVR模型预测结果与验证集实测相关性及RMSE分析(图3)

     test_simu = predict(model_init,
                          test)
      pl_point_df2 = data.frame(SIMU = as.numeric(test_simu),
                               TEST = test$shortwave)
      label_df2 = data.frame(
        x = c(-0.5,-0.5,-0.5),
        y = c(0.3,0.25,0.2),
        label = c(paste0("COR = ", round(cor(test_simu,test$shortwave),2)),
                  paste0("p = ", formatC(cor.test(test_simu,test$shortwave)$p.value)),
                  paste0("RMSE = ", round(rmse(test_simu,test$shortwave),3)))
      )
        ########
      p2_cor = ggplot()+
        geom_point(data = pl_point_df2,
                   aes(x = SIMU, y= TEST),
                   size = 5,
                   color = 'blue',
                   shape = 1)+
        geom_abline(intercept = 0,slope = 1,size = 1)+
        geom_text(data = label_df2,aes(x = x,y = y,label = label),size = 6,color = 'black',
                  hjust = 0)+
        theme_bw()+
        theme(
          axis.text =  element_text(face = 'bold',color = 'black',size = 12),
          axis.title =  element_text(face = 'bold',color = 'black',size = 14, hjust = .5)
        )+
        xlab('SIMULATION RESULT')+
        ylab('REAL SHORTWAVE')
      
      png('plot3.png',
          height = 20,
          width = 20,
          units = 'cm',
          res = 800)
      print(p2_cor)
      dev.off()
    
    图3 SVR模型预测结果与验证集实测相关性及RMSE分析

    8. SVR模型预测结果与验证集实测数据时间路径分析(图4)

    pl_line_df2 = data.frame(
        Date = date_test,
        SIMU = test_simu,
        REAL = test$shortwave
      )
      pl_line_df2 = melt(pl_line_df2,'Date')
      
      p1_line2 = ggplot()+
        geom_line(data = pl_line_df2,
                  aes(x = Date, y= value,color= variable,linetype = variable),
                  size= 1)+
        scale_color_manual(values = c('green','blue'))+
        scale_linetype_manual(values = c('solid','dashed'))+
        theme_bw()+
        theme(
          axis.text =  element_text(face = 'bold',color = 'black',size = 12),
          axis.title =  element_text(face = 'bold',color = 'black',size = 14, hjust = .5),
          legend.text = element_text(face = 'bold',color = 'black',size = 12),
          legend.title = element_blank(),
          legend.position = 'bottom',
          legend.direction = 'horizontal'
        )+
        xlab('Date')+
        ylab('Indices')
      
      png('plot4.png',
          height = 10,
          width = 20,
          units = 'cm',
          res = 800)
      print(p1_line2)
      dev.off()
    
    图4 SVR模型预测结果与验证集实测数据时间路径分析

    9. 模型模拟残差分析(图5)

    由图5可见,SVR模型在训练集与验证集的模拟残差大致服从于均值为0, 标准差为0.1 的正态分布,证明模型模拟与预测能力强。

    TRAIN_RES = model_init$fitted - train$shortwave
     TEST_RES = test_simu - test$shortwave
     
     df_res1 = data.frame(RES = TRAIN_RES, Type = "TRAIN_RES")
     df_res2 = data.frame(RES = TEST_RES, Type = "TEST_RES")
     df_res = rbind(df_res1,df_res2)
     label_df_res = data.frame(
       x = c(0.3,0.3,0.3,0.3),
       y = c(500,450,400,350),
       labels = c(
         paste0('TRAIN_MEAN = ',round(mean(TRAIN_RES),2)),
         paste0('TRAIN_SD = ',round(sd(TRAIN_RES),2)),
         paste0('TEST_MEAN = ',round(mean(TEST_RES),2)),
         paste0('TEST_SD = ',round(sd(TEST_RES),2))
       )
     )
     
     p5 = ggplot()+
       geom_density(data = df_res,aes(x = RES, y = ..count..,fill = Type))+
       geom_text(data = label_df_res,aes(x = x, y = y, label = labels),
                 color = 'black',fontface= 'bold',size = 5,hjust= 1)+
       theme_bw()+
       theme(
         axis.text =  element_text(face = 'bold',color = 'black',size = 12),
         axis.title =  element_text(face = 'bold',color = 'black',size = 14, hjust = .5),
         legend.text = element_text(face = 'bold',color = 'black',size = 12),
         legend.title = element_blank(),
         legend.position = 'bottom',
         legend.direction = 'horizontal'
       )+
       xlab('RES')+
       ylab('Count')
     
     png('plot5.png',
         height = 20,
         width = 20, 
         units = 'cm',
         res= 800)
     print(p5)
     dev.off()
    
    图5 模型残差分析

    10. 总结

    本篇主要解决了以下问题:
    1. 如何在R中构建核支持向量机模型?
    2. 如何对其预测结果进行定量化评估?

    11. 本篇所使用的软件包(没有的需要通过install.packages('包名')进行安装)

    library(reshape2)
    library(ggplot2)
    library(e1071)
    library(hydroGOF)
    

    12. 致谢

    首先,感谢大家的持续关注,小编会继续努力,持续更新下去的!

    大家如果觉得有用,还麻烦大家转发点赞加关注哈,也可以扩散到朋友圈,多谢大家啦~

    大家如果在使用本文代码的过程有遇到问题的,可以留言评论,也可以私信我哈~~

    最后,祝大家即将到来的元旦节快乐!

    小编联系方式

    相关文章

      网友评论

        本文标题:R - 核支持向量机预测模型构建详解:以短波辐射通量shortw

        本文链接:https://www.haomeiwen.com/subject/etsmoctx.html