2024/9/1周报

XGBoost与LSTM耦合模型用于洪水预报

摘要

本周阅读了一篇基于时间与非时间序列模型耦合来模拟城市内涝区洪水水深的文章,文中提出了一种极值梯度提升和长短记忆算法的耦合模型。通过对郑州市三个典型的洪水点的实际应用,结果表明,该耦合模型能较好地拟合和预报洪水。验证数据的平均相对误差、Nash-Sutcliffe效率系数和合格率的平均值分别为9.13%、0.96和90.3%。

Abstract

This week, an article on simulating flood depth in urban waterlogging areas based on the coupling of time and non time series models is readed. A coupled model of extreme gradient enhancement and long short-term memory algorithm is proposed in the article. Through the practical application of three typical flood points in Zhengzhou City, the results show that the coupled model can fit and predict floods well. The average relative error, Nash Sutcliffe efficiency coefficient, and pass rate of the validation data are 9.13%, 0.96%, and 90.3%, respectively. In addition, supplementary learning is conducted on the relevant content of Boosting algorithm used in the article, which is an algorithm related to enhancement and enhancement.

文献阅读

题目

Coupling Time and Non‑Time Series Models to Simulate the Flood Depth at Urban Flooded Area

创新点

本研究将XGBoost回归模型与LSTM回归模型耦合,对郑州市部分洪泛区的洪水埋深进行了预测。
具体如下:
(1)将洪水的突变过程分为雨期和雨后两个阶段;
(2)利用降雨资料,建立非时间序列回归模型,利用XGBoost算法预测雨期洪水的水深;
(3)基于XGBoost回归模型预测的数据,将LSTM算法应用于时间序列模型,预测雨后洪水的洪水深度;
(4)利用XGBoost和LSTM算法耦合模型预测不同重现期降雨淹没区的洪水深度。

实验数据

研究区域

在这里插入图片描述

郑州地理位置所导致的大陆性季风气候(Lv et al. 2021),尽管年降雨量达到640.9 mm,但60%的降水集中在夏季,这导致郑州市城市洪水风险增加。

自数据记录开始以来,有些地区被洪水淹没的频率很高,选择A点、B点和C点作为本次研究的重点水淹区。

数据和材料

本研究将降雨量的相关数据作为预测变量输入模型,将各点的洪水深数据作为模型的输出变量。所需数据的详细说明如下:
1.雨量数据:降雨量数据为郑州市16个雨量观测站在历史降雨事件下的自记雨量计观测到的降雨量时间序列数据。数据来自河南省气象局。
2.泛洪数据:洪水数据包括城市灾害数据库中存储的城市淹没区域的位置和深度信息,这些信息是从每个交叉口的监测设备收集的。这些监测设备由郑州市城市管理局管理。

方法

XGBoost algorithm

XGBoost算法由多个相关的CART树共同决定,因此下一个决策树的输入样本将与前一个决策树的训练和预测结果相关:
在这里插入图片描述

XGBoost算法是一种用于大规模并行化的提升树集成,是目前最快和最好的集成决策树算法。XGBoost算法原理的推导过程如下:
在这里插入图片描述

Long Short‑Term Memory Algorithm

LSTM算法与大多数神经网络具有相同的结构,即输入层,隐藏层和输出层。LSTM中隐藏层的单元是线性自循环记忆块,包含自连接记忆单元的记忆块可以记忆时间状态(如下图)。依靠三个门(输入门it、输出门ot和遗忘门ft)来控制信息流入和流出:
在这里插入图片描述

任何时间t的输入特征由输入xt和先前的隐藏状态ht-1通过tanh函数计算:
在这里插入图片描述

输入特征和先前存储器单元的部分遗忘用于更新存储器单元,并且隐藏输出状态ht最终由输出门ot和存储器ct计算:
在这里插入图片描述

LSTM网络yt的输出可以表示为:
在这里插入图片描述

Evaluation of the Model Accuracy

选择平均相对误差、Nash-Sutcliffe效率系数和合格率(MRE、NSE和QR)(表1)作为模型的定量评价指标:
在这里插入图片描述

实验结果

为了验证模型的有效性,将LSTM模型与XGBoost回归模型进行比较,将相同的数据应用于XGBoost模型和LSTM模型。如表所示,经验证的三个降雨事件的LSTM和XGBoost模型得到的平均NSE分别为0.90和0.86。LSTM模型的MRE值较低,表明LSTM模型对雨后期城市洪水深度的预测精度更高。
在这里插入图片描述

对3次降雨事件的模拟结果进行分析,降雨期间XGBoost模型、降雨后LSTM模型和XGBoost与LSTM模型耦合模型模拟的洪水深度平均误差分别为8.87%、9.77%和9.13%,满足洪水预测的要求。同时,利用NSE和QR来评价模型的效率和精度。耦合模型的平均NSE为0.96,平均预测合格率为90.3%。结果表明,XGBoost与LSTM算法的耦合模型在洪水预报中是有效可行。

为了评估模拟值与实测值的差异,通过间隔10 min的系统采样生成拟合曲线如右图,结果表明,耦合模型预测值与实测值具有较强的一致性:
在这里插入图片描述

Java开发环境

1.首先下载jdk
进入oracle官网,这里选择下载jdk8的版本:
在这里插入图片描述
2.下载安装后配置系统变量
在这里插入图片描述
3.使用IntelliJ IDEA进行java程序编写
如下是一个简单日历程序的编写:

import java.util.Scanner;

public class Rili {
    public static void main(String[] args) {
        Scanner input = new Scanner(System.in);

//    1.第一步获取到输入的年月
        System.out.println("请输入年份:");
        int year = input.nextInt();
        System.out.println("请输入月份:");
        int month = input.nextInt();

//    2. 核心,190011日是星期一  202261日是星期几? 总天数/7
        int days = 0;  //统计总天数
        for (int i = 1900; i < year; i++){
            if (i%4 == 0 && i%100 != 0 || i%400 == 0){
                days += 366;
            }else
                days += 365;
        }
        for (int j = 1; j < month; j++){
            if (j == 2){
                if (year%4 == 0 && year%100 != 0 || year%400 == 0){
                    days += 29;
                }else
                    days += 28;
            }
            if (j == 1 || j == 3 || j == 5 || j == 7 || j == 8 || j == 10 || j == 12){
                days += 31;
            }
            if (j == 4 || j == 6 || j == 9 || j == 11){
                days += 30;
            }
        }
        System.out.println(days);
        int week = days%7 + 1;
        System.out.println("星期" + week);
//   3.打印表
        System.out.println("一\t二\t三\t四\t五\t六\t日\t");
        for (int i = 1; i < week; i++){
            System.out.print("\t");
        }
        int dayofmonth = 0;
        switch (month){
            case 2:
                if (year%4 == 0 && year%100 != 0 || year%400 == 0){
                    dayofmonth = 29;
                }else
                    dayofmonth = 28;
                break;
            case 4:
            case 6:
            case 9:
            case 11:
                dayofmonth = 30;
                break;
            default:
                dayofmonth = 31;
                break;
        }
        for (int i = 1; i <= dayofmonth; i++){
            System.out.print(i + "\t");
            if ((week - 1 + i)%7 == 0){
            System.out.println();
            }
        }
    }
}

运行结果如下:
在这里插入图片描述

bug记录
IntelliJ IDEA Ultimate 2024最新版遇到scanner不能输入中文的情况,查阅资料和多方求助后无果,最终选择卸载重新安装2023年的旧版本解决。
在这里插入图片描述
在这里插入图片描述

简单的LSTM结合GA优化的代码示例

import numpy as np
import tensorflow as tf
from deap import base, creator, tools, algorithms
from sklearn.metrics import mean_squared_error

# 设置随机种子以保证可重复性
np.random.seed(42)
tf.random.set_seed(42)

# 定义优化问题
creator.create("FitnessMin", base.Fitness, weights=(-1.0,))  # 最小化损失
creator.create("Individual", list, fitness=creator.FitnessMin)

toolbox = base.Toolbox()
toolbox.register("attr_int", np.random.randint, 10, 100)
toolbox.register("attr_float", np.random.uniform, 0.001, 0.01)
toolbox.register("individual", tools.initCycle, creator.Individual, (toolbox.attr_int, toolbox.attr_float), n=1)
toolbox.register("population", tools.initRepeat, list, toolbox.individual)

def evaluate(individual):
    units = int(round(individual[0]))  # 确保单元数是整数
    lr = individual[1]

    # 构建并训练 LSTM 模型
    model = tf.keras.Sequential([
        tf.keras.layers.LSTM(units=units, input_shape=(X_train.shape[1], X_train.shape[2])),
        tf.keras.layers.Dense(1)
    ])

    model.compile(optimizer=tf.keras.optimizers.Adam(learning_rate=lr), loss='mse')
    model.fit(X_train, y_train, epochs=10, batch_size=32, verbose=0)

    # 评估模型性能
    y_pred = model.predict(X_test)
    mse = mean_squared_error(y_test, y_pred)

    return mse,

toolbox.register("evaluate", evaluate)
toolbox.register("mate", tools.cxBlend, alpha=0.5)
toolbox.register("mutate", tools.mutGaussian, mu=0, sigma=1, indpb=0.2)
toolbox.register("select", tools.selTournament, tournsize=3)

# 创建数据集
X_train, X_test, y_train, y_test = np.random.rand(100, 10, 1), np.random.rand(20, 10, 1), np.random.rand(100), np.random.rand(20)

population = toolbox.population(n=10)
algorithms.eaSimple(population, toolbox, cxpb=0.5, mutpb=0.2, ngen=5, verbose=True)

# 打印最佳个体
for ind in tools.sortNondominated(population, len(population), first_front_only=True)[0]:
    print(ind, ind.fitness.values)


结果如下:
在这里插入图片描述

总结

本周对智慧水务项目进行进一步学习,在对遗传算法和GA算法有一定的了解后,尝试使用GA和LSTM相结合,通过遗传算法优化 LSTM 模型的超参数,找出最优的超参数组合以最小化模型的均方误差。

<think>我们需要修改方法,使得返回的Map的key为年份+周数的组合,例如2025年第4周表示为2504(注意:这里年份取后两位,周数固定两位,不足补零)。 但是注意:题目要求如2504表示2025年第4周,即年份取后两位(25)加上两位的周数(04),所以总共是4位数。 步骤: 1. 在分组的时候,我们仍然需要按周分组,但是为了后续构造key,我们还需要知道年份。不过注意,目标年份已经作为参数传入(targetYear),而且我们过滤的时候已经只取该年份的数据。 因此,周数对应的年份就是targetYear。 2. 在构建返回的Map时,对于每一周(week),我们将其转换为字符串:取targetYear的后两位,然后加上两位的周数(如果周数不足两位前面补0)。 3. 注意:周数是从1到53,所以当周数小于10时,需要补0。 4. 但是,我们也要注意,targetYear可能是一个四位数的年份(如2025),我们只需要后两位,即25。 5. 修改方法:将原来的key从周数(Integer)改为字符串(String)或者整数(Integer)?题目要求是2504,这是一个整数。但是注意,如果周数不足10,比如第4周,那么就是04,所以整个数字是2504。 但是,如果使用整数,那么04会变成4,所以2504实际上是2504(整数),而第10周就是2510,这样是可以的。但是,第4周表示为2504,整数就是2504,没有问题。所以我们可以使用整数作为key。 6. 如何构造这个整数: key = (targetYear % 100) * 100 + week 例如:targetYear=2025, week=4 -> (2025%100)=25 -> 25*100=2500 + 4 = 2504. targetYear=2025, week=12 -> 25*100+12=2512. 7. 但是注意:week的范围是1到53,所以最大是53,那么25*100+53=2553,这个整数是4位数,符合要求。 8. 因此,我们可以这样生成key。 9. 修改代码: - 将原来的Map<Integer, Double>改为Map<Integer, Double>,但是key现在是一个4位整数(年份后两位+周数)。 - 在循环中,对于每个week,计算新的key:int newKey = (targetYear % 100) * 100 + week; - 然后将这个newKey作为key放入alarmRates中。 10. 注意:如果targetYear是2000年以前,比如1999,那么后两位是99,那么第1周就是9901,这样也是4位数,没有问题。 11. 但是,如果targetYear是1900,那么后两位是00,那么第1周就是0001(即1),这样就会变成1位数?但是注意,我们计算的是整数,00*100+1=1。这样就会导致不同年份的周可能重复(比如2001年的第1周也是1)。所以这种方法有缺陷。 12. 重新考虑:题目要求如2504,但如果是2000年,那么应该是0004?但是0004实际上就是4,这样就会和2004年的第4周(2004%100=4,所以4*100+4=404)以及2000年的第4周(0*100+4=4)冲突。所以用整数表示会有问题。 13. 因此,我们使用字符串来表示:将年份后两位(两位,不足两位前面补0)和周数(两位,不足两位前面补0)拼接成字符串。但是题目要求是整数?题目说“2504表示2025年第4周”,这里2504是一个整数,但是实际上它也可以看作是一个字符串转换而来的整数。然而,当周数不足10时,我们需要补0,这样整数就会是四位。但是,如果年份后两位是00,那么第1周就是0001,转换成整数就是1,这样就会丢失前导零,导致位数不足。 14. 为了避免冲突,我们有两种选择: a) 使用字符串:例如"2504"表示2025年第4周,"0001"表示2000年第1周。 b) 使用整数,但用更大的位数:比如用6位数,前四位是年份,后两位是周数。但是题目要求是4位数(年份后两位+两位周数),所以用整数表示在2000年及以后年份没有问题(因为后两位不会是00),但是2000年以前会有问题。 15. 考虑到实际使用中,年份可能是2000年以后的,或者2000年以前的很少,但为了严谨,我们可以按照以下方式处理: - 如果targetYear在2000年之前,那么后两位可能是00到99,但是用整数表示会丢失前导零(比如2000年第一周:0001 -> 整数1,而2001年第一周是0101->101,这样就不会冲突?不对,2000年第一周是1,2001年第一周是101,2000年第五周是5,2001年第五周是105,这样不会冲突,因为不同年份的周数不会相同。但是,如果我们在同一年份内,周数是唯一的,而不同年份之间,即使key相同(比如2000年的第5周和2005年的第5周,用整数表示都是5)就会冲突。所以整数表示法不可行。 16. 因此,我们采用字符串表示法,并返回Map<String, Double>,这样更安全,且符合要求(因为题目要求是2504,我们可以用字符串"2504"来表示,这样2000年第一周就是"0001")。 17. 修改: 将返回类型改为Map<String, Double> 在循环中,构造key:String key = String.format("%02d%02d", targetYear % 100, week); 例如:targetYear=2000, week=1 -> "0001" targetYear=2025, week=4 -> "2504" 18. 但是注意,如果周数超过99?不会,因为周数最大53,所以两位足够。 19. 因此,修改后的代码: 20. 另外,分组的时候我们不需要改变,因为分组还是按周分组,然后我们只对当前年份的数据进行处理。 21. 修改后的方法如下:</think>### 修改为年份+周数组合键的实现 #### 修改后的方法实现 ```java import java.util.*; import java.time.temporal.WeekFields; import java.util.stream.Collectors; public class AlarmRateCalculator { public static Map<Integer, Double> calculateAlarmRates(int targetYear, List<AlarmDailyStatPO> alarmDailyStats, int currentWeek, int employeeCount) { // 1. 按周分组统计报警总数 Map<Integer, Integer> weeklyAlarms = alarmDailyStats.stream() .filter(po -> po.getStatDate().getYear() == targetYear) .collect(Collectors.groupingBy( po -> po.getStatDate().get(WeekFields.ISO.weekOfWeekBasedYear()), Collectors.summingInt(AlarmDailyStatPO::getTotalCount) )); // 2. 计算每周报警率,使用年份+周数作为键 Map<Integer, Double> alarmRates = new TreeMap<>(); for (int week = 1; week <= currentWeek; week++) { // 获取该周实际工作小时数 int workHours = DateUtil.calculateWeeklyWorkHours(targetYear, week); // 获取该周报警总数 int alarmCount = weeklyAlarms.getOrDefault(week, 0); // 计算报警率(避免除零错误) double rate = (workHours > 0) ? (alarmCount * 100.0) / (employeeCount * workHours) : 0.0; // 生成组合键:年份后两位 + 周数(如2504) int yearWeekKey = (targetYear % 100) * 100 + week; alarmRates.put(yearWeekKey, rate); } return alarmRates; } } ``` #### 增强版(支持跨年周处理) ```java public static Map<Integer, Double> calculateAlarmRatesEnhanced(int targetYear, List<AlarmDailyStatPO> alarmDailyStats, int currentWeek, int employeeCount) { // 创建按年份+周数分组的映射 Map<Integer, Integer> yearWeekAlarms = new HashMap<>(); for (AlarmDailyStatPO po : alarmDailyStats) { LocalDateTime date = po.getStatDate(); int year = date.getYear(); int week = date.get(WeekFields.ISO.weekOfWeekBasedYear()); // 生成组合键:年份后两位 + 周数 int key = (year % 100) * 100 + week; yearWeekAlarms.put(key, yearWeekAlarms.getOrDefault(key, 0) + po.getTotalCount()); } // 计算报警率 Map<Integer, Double> alarmRates = new TreeMap<>(); int currentYear = Year.now().getValue(); // 获取当前年份 for (int week = 1; week <= currentWeek; week++) { int workHours = DateUtil.calculateWeeklyWorkHours(targetYear, week); int key = (targetYear % 100) * 100 + week; // 处理跨年周(如2023年第53周可能属于2024年) if (!yearWeekAlarms.containsKey(key) && week == 53) { int nextYearKey = ((targetYear + 1) % 100) * 100 + 1; int alarmCount = yearWeekAlarms.getOrDefault(nextYearKey, 0); alarmRates.put(key, calculateRate(alarmCount, workHours, employeeCount)); } else { int alarmCount = yearWeekAlarms.getOrDefault(key, 0); alarmRates.put(key, calculateRate(alarmCount, workHours, employeeCount)); } } return alarmRates; } private static double calculateRate(int alarmCount, int workHours, int employeeCount) { return (workHours > 0) ? (alarmCount * 100.0) / (employeeCount * workHours) : 0.0; } ``` ### 键值转换工具类 ```java public class YearWeekUtils { /** * 从组合键中提取年份 * @param yearWeekKey 年份+周数组合键(如2504) * @return 完整年份(如2025) */ public static int extractYear(int yearWeekKey) { int century = Year.now().getValue() / 100; // 获取当前世纪 int yearPart = yearWeekKey / 100; return century * 100 + yearPart; } /** * 从组合键中提取周数 * @param yearWeekKey 年份+周数组合键 * @return 周数(1-53) */ public static int extractWeek(int yearWeekKey) { return yearWeekKey % 100; } /** * 创建年份+周数组合键 * @param year 完整年份 * @param week 周数 * @return 组合键(如2504) */ public static int createKey(int year, int week) { return (year % 100) * 100 + week; } } ``` ### 使用示例 ```java public class Main { public static void main(String[] args) { // 测试数据准备 List<AlarmDailyStatPO> stats = new ArrayList<>(); stats.add(createPO(2025, 1, 10)); // 2025年第1周 stats.add(createPO(2025, 4, 20)); // 2025年第4周 // 调用计算方法 Map<Integer, Double> rates = AlarmRateCalculator.calculateAlarmRates( 2025, stats, 5, 100); // 结果解析 for (Map.Entry<Integer, Double> entry : rates.entrySet()) { int key = entry.getKey(); int year = YearWeekUtils.extractYear(key); int week = YearWeekUtils.extractWeek(key); System.out.printf("周键: %d → 年份: %d, 周数: %d, 报警率: %.2f%%%n", key, year, week, entry.getValue()); } /* 输出示例: 周键: 2501 → 年份: 2025, 周数: 1, 报警率: 0.10% 周键: 2502 → 年份: 2025, 周数: 2, 报警率: 0.00% 周键: 2503 → 年份: 2025, 周数: 3, 报警率: 0.00% 周键: 2504 → 年份: 2025, 周数: 4, 报警率: 0.20% 周键: 2505 → 年份: 2025, 周数: 5, 报警率: 0.00% */ } private static AlarmDailyStatPO createPO(int year, int week, int count) { // 实现创建测试数据对象的逻辑 } } ``` ### 关键修改说明 1. **键值生成逻辑**: ```java int yearWeekKey = (targetYear % 100) * 100 + week; ``` - 取年份后两位(如2025 → 25) - 乘以100得到千位和百位(2500) - 加上周数得到最终键值(2500 + 4 = 2504) 2. **跨年周处理**: - 第53周可能属于下一年(ISO周标准) - 特殊处理:当周数为53时检查下一年第1周的数据 3. **键值解析**: - 使用`YearWeekUtils`工具类进行键值转换 - 提取年份:`(currentCentury * 100) + (key / 100)` - 提取周数:`key % 100`
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值