使用java实现蒙特卡洛模拟风险预测功能

巧合下看到了蒙特卡洛模拟风险预测，比较感兴趣，这里使用SpringBoot + Vue3来实现这个功能。

注意因为本人没使用过水晶球这款软件，所以这里的算法处理以及功能也是基于博客，以及网上的一些资料做的，如果有疏漏或者是不对的地方，可以在评论区说明，有时间的话我再改。

参照的文章：

Crystal Ball—甲骨文水晶球风险管理软件（概念以及实战——基础案例篇）_crystal ball软件-CSDN博客

https://zhuanlan.zhihu.com/p/477816011

项目代码是基于开源项目 小诺/Snowy 做的，全部技术栈 SpringBoot+MybatisPlus+ VUE3。

当前项目代码放在 李浩/snowy-interesting

蒙特卡洛模拟

蒙特卡罗（Monte Carlo）方法，又称随机抽样或统计试验方法。基本思想是通过某种“试验”的方法，得到事件出现的频率，或者这个随机变数的平均值，并用它们作为问题的的解。蒙特卡罗解题归结为三个主要步骤：构造或描述概率过程；实现从已知概率分布抽样；建立各种估计量。

说实话刚看的时候实在有点懵逼，尤其是各种计算公式头皮发麻，看了下各种博客后，自己的理解是现在需要获得某件事的成功的概率，事件中包含多个值，根据每个值获得随机值，然后按照指定的运算获得最终值，查看最终值是否在预料的结果内，如果在则说明成功，否则记作失败，循环数次这个操作，其中成功所占比就是本次蒙特卡洛模拟的结果（这里只是自己的理解，错的话可以在评论区指出）。

简单点讲，现在有事件y，y事件中有a,b两个值，y = a+b，a的话大于0小于100，b大于20小于40，然后根据a和b的取值范围生成一个随机值，然后拿着随机值求y，如果y的值小于100则表示这个结果正确，大于100则表示计算错误，然后循环10000次，最后10000个结果中为正确的数量所占比，就是本次模拟的结果，当然也可以有更加复杂的运算处理。

举个现实生活例子：

现在 要计算整个项目工时是否会超时完成的概率。整个活动的完成时间受设计，开发，测试三个值影响，其中设计步骤是最大值20，最小是8的正态分布值，开发是最大值是40最小为10的正态分布值，测试是最大值是20最小值5的正态分布值；

根据设计，开发，测试的分布开始获得三个值，而活动最终时间 = 设计+开发+测试+12。根据公式循环10000次然后计算出10000个结果，结果在平均值以上并且小于100为正常完成，通过运算可以获得1000次成功，也就是项目正常完成的概率是1000 / 10000 = 0.1 = 10%的概率，那么也就意味我我们需要更改每个阶段或者值的范围来达到理想结果。

从这两个例子中其实就可以理解蒙特卡罗解题步骤：

1、构造或描述概率过程 -》 确定最终结果以及参数的计算方式

2、实现从已知概率分布抽样 -》 按照分布抽取参数的随机值

3、建立各种估计量 -》 运算，然后获得概率，然后重复1,2直到理想结果

当然蒙特卡罗模拟不过只是如此简单，还有更加复杂的处理方式，这里不作为讨论

Python实现的蒙特卡洛模拟

这里是使用Python来实现蒙特卡洛模拟，这里的代码主要用于理解蒙特卡洛模拟，而不是功能模块，代码主要含义是假设现在有一个圆，然后随机生成位置，模拟10000，在圆内则表示击中

# 设置飞镖数量
num_darts=10000
# 击中数
darts_in_circle=0

# 存储飞镖坐标
x_coords_in, y_coords_in, x_coords_out, y_coords_out= [], [], [], []

# 设置图形数量
num_figures=6
darts_per_figure=num_darts//num_figures

# 生成单位圆
theta=np.linspace(0, 2*np.pi, 100)
unit_circle_x=np.cos(theta)
unit_circle_y=np.sin(theta)

pi_estimates = []

# 模拟投掷飞镖
for i in range(num_darts):
    # 获得随机值
    x, y=random.uniform(-1, 1), random.uniform(-1, 1)

    # 计算非负数的平方根
    # ** 幂运算
    if math.sqrt(x**2+y**2) <=1:
        darts_in_circle+=1
        x_coords_in.append(x)
        y_coords_in.append(y)
    else:
        x_coords_out.append(x)
        y_coords_out.append(y)

    if (i+1) % darts_per_figure == 0:
        pi_estimate=4*darts_in_circle/ (i+1)

        pi_estimates.append(pi_estimate)
        
        estimated_circle_radius= pi_estimate / 4
        estimated_circle_x= estimated_circle_radius * np.cos(theta)
        estimated_circle_y= estimated_circle_radius * np.sin(theta)

        fig=go.Figure()


        fig.add_trace(
            go.Scattergl(
                x=x_coords_in, 
                y=y_coords_in, 
                mode='markers', 
                name='Darts Inside Circle', 
                marker=dict(color='green', size=4, opacity=0.8))
        )
        fig.add_trace(
            go.Scattergl(
                x=x_coords_out, 
                y=y_coords_out, 
                mode='markers', 
                name='Darts Outside Circle', 
                marker=dict(color='red', size=4, opacity=0.8))
        )


        fig.add_trace(go.Scatter(x=unit_circle_x, y=unit_circle_y, mode='lines', name='Unit Circle', line=dict(color='green', width=3)))
        fig.add_trace(go.Scatter(x=estimated_circle_x, y=estimated_circle_y, mode='lines', name='Estimated Circle', line=dict(color='blue', width=3)))


        fig.update_layout(title=f"Figure {chr(97+ (i+1) //darts_per_figure-1)}: Thrown Darts: {(i+1)}, Estimated Pi: {pi_estimate}", width=600, height=600, xaxis=dict(constrain="domain", range=[-1, 1]), yaxis=dict(scaleanchor="x", scaleratio=1, range=[-1, 1]), legend=dict(yanchor="top", y=0.99, xanchor="left", x=0.01))

        fig.show()

        pio.write_image(fig, f"fig2{chr(97+ (i+1) //darts_per_figure-1)}.png")

最后的运行结果，其中绿点标识击中

java实现的蒙特卡洛模拟

代码基于snowy开源项目做的，这里的项目不是简单的模拟，而是功能模块，前端代码在snowy-admin-web项目中monte文件夹下，后端在interesting-plugin-monte项目中

1、数据库设计

monte_simulate 模拟表，存放模拟相关数据

CREATE TABLE `monte_simulate`  (
  `ID` varchar(20) CHARACTER SET utf8 COLLATE utf8_general_ci NOT NULL COMMENT 'ID',
  `NAME` varchar(255) CHARACTER SET utf8 COLLATE utf8_general_ci NULL DEFAULT NULL COMMENT '名称',
  `RESULT_MAX` double(255, 0) NULL DEFAULT NULL COMMENT '结果的最大值',
  `RESULT_MIN` double(255, 0) NULL DEFAULT NULL COMMENT '结果的最小值',
  `CORRECT_ODDS` varchar(255) CHARACTER SET utf8 COLLATE utf8_general_ci NULL DEFAULT NULL COMMENT '正确概率',
  `SIMULATE_NUM` int(11) NULL DEFAULT NULL COMMENT '模拟次数',
  `FORMULA` varchar(255) CHARACTER SET utf8 COLLATE utf8_general_ci NULL DEFAULT NULL COMMENT '计算公式',
  `CREATE_TIME` datetime(0) NULL DEFAULT NULL COMMENT '创建时间',
  `CREATE_USER` varchar(20) CHARACTER SET utf8mb4 COLLATE utf8mb4_general_ci NULL DEFAULT NULL COMMENT '创建用户',
  `UPDATE_TIME` datetime(0) NULL DEFAULT NULL COMMENT '修改时间',
  `UPDATE_USER` varchar(20) CHARACTER SET utf8mb4 COLLATE utf8mb4_general_ci NULL DEFAULT NULL COMMENT '修改用户',
  PRIMARY KEY (`ID`) USING BTREE
) ENGINE = InnoDB CHARACTER SET = utf8 COLLATE = utf8_general_ci COMMENT = '蒙特卡洛模拟 ' ROW_FORMAT = Dynamic;

monte_simulate_params 模拟参数表，主要是写入本次模拟所使用的参数

CREATE TABLE `monte_simulate_params`  (
  `ID` varbinary(20) NOT NULL COMMENT 'ID',
  `SIMULATE_ID` varchar(20) CHARACTER SET utf8 COLLATE utf8_general_ci NULL DEFAULT NULL COMMENT '模拟ID',
  `TYPE` char(1) CHARACTER SET utf8 COLLATE utf8_general_ci NULL DEFAULT NULL COMMENT '类型，1变量，2常量',
  `DISTRIBUTION_TYPE` char(1) CHARACTER SET utf8 COLLATE utf8_general_ci NULL DEFAULT NULL COMMENT '分布类型， 1正态分布',
  `MEAN_VALUE` double(255, 3) NULL DEFAULT NULL COMMENT '均值',
  `SD_VALUE` double(255, 3) NULL DEFAULT NULL COMMENT '平方差值',
  `PARAM_KEY` varchar(255) CHARACTER SET utf8 COLLATE utf8_general_ci NULL DEFAULT NULL COMMENT '参数key',
  PRIMARY KEY (`ID`) USING BTREE
) ENGINE = InnoDB CHARACTER SET = utf8 COLLATE = utf8_general_ci COMMENT = '蒙特卡洛模拟 参数' ROW_FORMAT = Dynamic;

monte_calculate_log 模拟记录表，用于存储对本次模拟

CREATE TABLE `monte_calculate_log`  (
  `ID` varchar(20) CHARACTER SET utf8 COLLATE utf8_general_ci NOT NULL COMMENT 'ID',
  `SIMULATE_ID` varchar(20) CHARACTER SET utf8 COLLATE utf8_general_ci NULL DEFAULT NULL COMMENT '模拟ID',
  `SIMULATE_DATA` text CHARACTER SET utf8 COLLATE utf8_general_ci NULL COMMENT '模拟数据',
  `PARAMS` mediumtext CHARACTER SET utf8 COLLATE utf8_general_ci NULL COMMENT '相关参数',
  `RESULT` mediumtext CHARACTER SET utf8 COLLATE utf8_general_ci NULL COMMENT '计算结果',
  `TRUE_NUM` int(11) NULL DEFAULT NULL COMMENT '正确数量',
  `SIMULATE_NUM` int(11) NULL DEFAULT NULL COMMENT '模拟数量',
  `CREATE_TIME` datetime(0) NULL DEFAULT NULL COMMENT '创建时间',
  PRIMARY KEY (`ID`) USING BTREE
) ENGINE = InnoDB CHARACTER SET = utf8 COLLATE = utf8_general_ci COMMENT = '蒙特卡洛模拟 计算记录' ROW_FORMAT = Dynamic;

2、模块设计

设想是首先用户先创建模拟数据，定义本次模拟次数，正确值范围，计算公式以及名称等常规数据，然后定义模拟中使用的参数，定义参数类型分布方式（当前只写了正态分布的，其他的以后有时间再说），定义完成后开始模拟，每次模拟生成一条日志记录，记录中保存当前模拟的参数以及结果，然后基于结果使用echarts折线图显示。

基于大致思想，后端代码主要分为三部分

simulate 模拟

simulateParams 模拟参数

calculateLog 日志

每部分包含自己的controller，service，Mapper

3、主要代码

个人认为主要代码的话后端是读取计算公式，然后根据计算公式以及参数生成并计算

Integer trueNum = 0; // 正确数量
for(int a=0; a< simulateNum; a++){
    Map<String, Double> variables = new HashMap<>();
    for (MonteSimulateParams params : param.getParams()) {
        double randomValue = 0.0;
        if(params.getDistributionType().equals(SimulateParamsTypeEnum.VARIABLE.getValue())){
            // 根据正态分布生成随机值
            NormalDistribution distribution = new NormalDistribution(params.getMeanValue(),
                                                                     params.getSdValue());
            randomValue = distribution.sample();
        }else{
            // 常量
            randomValue = params.getMeanValue();
        }
        variables.put(params.getParamKey(), randomValue);
    }
    vStr.append(JSONUtil.toJsonStr(variables) + ",");
    double result = evaluate(monteSimulate.getFormula(), variables);
    rStr.append(result);
    rStr.append(",");

    if(result >= param.getResultMin() && result <= param.getResultMax()){
        trueNum++;
    }
}

解析计算公式

public double evaluate(String expression, Map<String, Double> variables) {
        // 去除所有空格
        String expr = expression.replaceAll("\\s+", "");

        // 分割表达式为tokens（变量名和运算符）
        List<String> tokens = splitTokens(expr);

        // 校验tokens的数量是否为奇数（例如，操作数、运算符、操作数...）
        if (tokens.size() % 2 == 0) {
            throw new IllegalArgumentException("无效的表达式格式");
        }

        // 将变量名替换为对应的值，并转换为操作数和运算符的列表
        List<Object> processedTokens = new ArrayList<>();
        for (String token : tokens) {
            if (isOperator(token)) {
                processedTokens.add(token);
            } else {
                Double value = variables.get(token);
                if (value == null) {
                    throw new IllegalArgumentException("变量 '" + token + "' 不存在");
                }
                processedTokens.add(value);
            }
        }

        // 处理乘除运算
        List<Object> afterMulDiv = processMulDiv(processedTokens);

        // 处理加减运算
        return processAddSub(afterMulDiv);
    }

// 使用正则表达式分割tokens
private static List<String> splitTokens(String expr) {
    List<String> tokens = new ArrayList<>();
    Pattern pattern = Pattern.compile("([a-zA-Z0-9]+)|([-+*/])");
    Matcher matcher = pattern.matcher(expr);
    while (matcher.find()) {
        String token = matcher.group();
        tokens.add(token);
    }
    return tokens;
}

// 判断是否为运算符
private static boolean isOperator(String token) {
    return token.matches("[-+*/]");
}

// 处理乘除运算
private static List<Object> processMulDiv(List<Object> tokens) {
    List<Object> result = new ArrayList<>();
    int i = 0;
    while (i < tokens.size()) {
        Object current = tokens.get(i);
        if (current instanceof String) {
            String op = (String) current;
            if (op.equals("*") || op.equals("/")) {
                // 取出前一个操作数
                if (result.isEmpty()) {
                    throw new IllegalArgumentException("无效的表达式：运算符在开头");
                }
                double left = (Double) result.remove(result.size() - 1);
                i++; // 移动到下一个元素（操作数）
                if (i >= tokens.size()) {
                    throw new IllegalArgumentException("无效的表达式：运算符后无操作数");
                }
                Object rightObj = tokens.get(i);
                if (!(rightObj instanceof Double)) {
                    throw new IllegalArgumentException("无效的表达式：运算符后不是数字");
                }
                double right = (Double) rightObj;
                double calc = calculate(left, right, op);
                result.add(calc);
                i++; // 处理完毕，i前进
            } else {
                // 加减运算符直接加入结果列表
                result.add(op);
                i++;
            }
        } else {
            result.add(current);
            i++;
        }
    }
    return result;
}

// 处理加减运算
private static double processAddSub(List<Object> tokens) {
    if (tokens.isEmpty()) {
        return 0.0;
    }
    if (tokens.size() % 2 == 0) {
        throw new IllegalArgumentException("无效的表达式：操作数与运算符不匹配");
    }
    double res = (Double) tokens.get(0);
    for (int i = 1; i < tokens.size(); i += 2) {
        String op = (String) tokens.get(i);
        double num = (Double) tokens.get(i + 1);
        res = calculate(res, num, op);
    }
    return res;
}

// 执行基础运算
private static double calculate(double a, double b, String op) {
    switch (op) {
        case "+":
            return a + b;
        case "-":
            return a - b;
        case "*":
            return a * b;
        case "/":
            if (b == 0) {
                throw new ArithmeticException("除以零错误");
            }
            return a / b;
        default:
            throw new IllegalArgumentException("未知运算符: " + op);
    }
}

前端主要有点意思的是echarts，之前都是直接用ID创建echarts对象，第一次在Vue3下弹窗里面使用，发现要使用ref否则会报错，这里记录下也

<a-form-item-rest>
  <div ref="basicLineChartRef" style="width: 100%; height: 500px"></div>
</a-form-item-rest>

const basicLineChartRef = ref()

const onOpen = (logId) => {
		open.value = true
		if (logId) {
      ...
        let Echarts = echarts.init(basicLineChartRef.value)
        chartOption.value = {
					tooltip: {}, // 默认提示框配置
					xAxis: {
						data: xDataArr // X 轴数据
					},
					yAxis: {}, // Y 轴默认配置
					series: [
						{
							type: 'line', // 图表类型为柱状图
							data: data.simulateResult  // 数据值
						}
					]
				}
        Echarts.setOption(chartOption.value)
      ...
    }
}

代码也就写这么多吧，感兴趣的话可以去git仓库里看看，后期如果有哪些比较有意思的想法功能，也会后续写入添加。

最后再重申下，这个功能也是偶然了解到后，对这个功能比较感兴趣，想实际试试实际项目功能是要怎么处理的，而日常工作中也并没有用过水晶球工具以及相关功能，只是根据一些博客上了解到的流程处理，如果有不对的，可以在留言里面说明，后面再改。