Browse Source

DEA模型Example

environment_dev
gjh 4 weeks ago
parent
commit
f9e22734e0
  1. 64
      ruoyi-modules/ruoyi-demo/src/main/java/org/dromara/demo/controller/test/ExtendedDEAExample.java
  2. 551
      ruoyi-modules/ruoyi-demo/src/main/java/org/dromara/demo/controller/test/ExtendedDEAUtils.java

64
ruoyi-modules/ruoyi-demo/src/main/java/org/dromara/demo/controller/test/ExtendedDEAExample.java

@ -0,0 +1,64 @@
package org.dromara.demo.controller.test;
import java.util.List;
/**
* 扩展DEA模型BCC和非期望SBM使用示例
*/
public class ExtendedDEAExample {
public static void main(String[] args) {
// 示例1:BCC模型使用示例
System.out.println("===== BCC模型示例 =====");
// 3个DMU,2个投入指标,1个产出指标
double[][] inputsBCC = {
{2,2,1,2,1,1}, // DMU 0的投入
{3,2,5,4,4,1}, // DMU 0的投入
{2,3,3,3,5,3}, // DMU 0的投入
{5,4,1,2,2,1}, // DMU 0的投入
{2,3,3,3,5,3}, // DMU 0的投入
{5,4,1,2,2,1}, // DMU 0的投入
};
double[][] outputsBCC = {
{4}, // DMU 0的产出
{3}, // DMU 0的产出
{2}, // DMU 0的产出
{3}, // DMU 0的产出
{5}, // DMU 0的产出
{5}, // DMU 0的产出
};
List<ExtendedDEAUtils.DEAResult> bccResults =
ExtendedDEAUtils.calculateAllBCC(inputsBCC, outputsBCC);
for (ExtendedDEAUtils.DEAResult result : bccResults) {
System.out.println(result);
}
// 示例2:非期望SBM模型使用示例
System.out.println("\n===== 非期望SBM模型示例 =====");
// 3个企业(DMU),2个投入(资本、劳动),1个期望产出(利润),1个非期望产出(污染)
double[][] inputsSBM = {
{100, 50}, // 企业1投入:资本100,劳动50
{150, 60}, // 企业2投入
{80, 40} // 企业3投入
};
double[][] desirableOutputs = {
{200}, // 企业1利润
{240}, // 企业2利润
{150} // 企业3利润
};
double[][] undesirableOutputs = {
{30}, // 企业1污染排放
{40}, // 企业2污染排放
{20} // 企业3污染排放
};
List<ExtendedDEAUtils.DEAResult> sbmResults =
ExtendedDEAUtils.calculateAllUndesirableSBM(inputsSBM, desirableOutputs, undesirableOutputs);
for (ExtendedDEAUtils.DEAResult result : sbmResults) {
System.out.println(result);
}
}
}

551
ruoyi-modules/ruoyi-demo/src/main/java/org/dromara/demo/controller/test/ExtendedDEAUtils.java

@ -0,0 +1,551 @@
package org.dromara.demo.controller.test;
import org.apache.commons.math3.analysis.MultivariateFunction;
import org.apache.commons.math3.optim.*;
import org.apache.commons.math3.optim.linear.*;
import org.apache.commons.math3.optim.nonlinear.scalar.GoalType;
import org.apache.commons.math3.optim.nonlinear.scalar.MultivariateOptimizer;
import org.apache.commons.math3.optim.nonlinear.scalar.ObjectiveFunction;
import org.apache.commons.math3.optim.nonlinear.scalar.noderiv.NelderMeadSimplex;
import org.apache.commons.math3.optim.nonlinear.scalar.noderiv.SimplexOptimizer;
import java.util.ArrayList;
import java.util.List;
/**
* 扩展的DEA工具类包含BCC模型和非期望SBM模型
*/
public class ExtendedDEAUtils {
// ===================== BCC模型(规模报酬可变) =====================
/**
* 计算指定DMU的BCC模型效率值面向投入
* BCC模型在CCR基础上增加了sum(λ_j) = 1的约束考虑规模报酬可变
*/
public static DEAResult calculateBCC(double[][] inputs, double[][] outputs, int dmuIndex) {
// 验证输入数据合法性
validateInput(inputs, outputs, dmuIndex);
int n = inputs.length; // DMU数量
int m = inputs[0].length; // 投入指标数量
int s = outputs[0].length; // 产出指标数量
// 变量数量: n个lambda + 1个theta + m个投入松弛变量 + s个产出松弛变量
int varCount = n + 1 + m + s;
// 目标函数: 最小化theta(第n个变量)
double[] objectiveCoefficients = new double[varCount];
objectiveCoefficients[n] = 1.0; // theta的系数为1,其他变量系数为0
// 创建线性规划问题
LinearObjectiveFunction objective = new LinearObjectiveFunction(objectiveCoefficients, 0);
// 约束条件列表
List<LinearConstraint> constraints = new ArrayList<>();
// 1. 投入约束: sum(lambda_j * x_ij) + s_i^- = theta * x_ik (i=1..m)
for (int i = 0; i < m; i++) {
double[] coefficients = new double[varCount];
// 设置lambda_j的系数
for (int j = 0; j < n; j++) {
coefficients[j] = inputs[j][i];
}
// 设置theta的系数
coefficients[n] = -inputs[dmuIndex][i];
// 设置投入松弛变量的系数
coefficients[n + 1 + i] = 1.0;
constraints.add(new LinearConstraint(coefficients, Relationship.EQ, 0.0));
}
// 2. 产出约束: sum(lambda_j * y_rj) - s_r^+ = y_rk (r=1..s)
for (int r = 0; r < s; r++) {
double[] coefficients = new double[varCount];
for (int j = 0; j < n; j++) {
coefficients[j] = outputs[j][r];
}
coefficients[n + 1 + m + r] = -1.0;
constraints.add(new LinearConstraint(coefficients, Relationship.EQ, outputs[dmuIndex][r]));
}
// 3. 变量非负约束
for (int j = 0; j < n; j++) {
double[] coefficients = new double[varCount];
coefficients[j] = 1.0;
constraints.add(new LinearConstraint(coefficients, Relationship.GEQ, 0.0));
}
// 投入松弛变量非负
for (int i = 0; i < m; i++) {
double[] coefficients = new double[varCount];
coefficients[n + 1 + i] = 1.0;
constraints.add(new LinearConstraint(coefficients, Relationship.GEQ, 0.0));
}
// 产出松弛变量非负
for (int r = 0; r < s; r++) {
double[] coefficients = new double[varCount];
coefficients[n + 1 + m + r] = 1.0;
constraints.add(new LinearConstraint(coefficients, Relationship.GEQ, 0.0));
}
// 4. theta <= 1约束
double[] thetaCoeff = new double[varCount];
thetaCoeff[n] = 1.0;
constraints.add(new LinearConstraint(thetaCoeff, Relationship.LEQ, 1.0));
// 5. BCC模型特有约束:sum(lambda_j) = 1(规模报酬可变)
double[] sumLambdaCoeff = new double[varCount];
for (int j = 0; j < n; j++) {
sumLambdaCoeff[j] = 1.0;
}
constraints.add(new LinearConstraint(sumLambdaCoeff, Relationship.EQ, 1.0));
// 求解线性规划问题
SimplexSolver solver = new SimplexSolver();
PointValuePair solution = solver.optimize(
new LinearConstraintSet(constraints),
objective,
GoalType.MINIMIZE,
new NonNegativeConstraint(false),
new MaxIter(1000)
);
// 解析结果
double[] solutionPoint = solution.getPoint();
double theta = solutionPoint[n]; // 效率值
// 投入松弛变量
double[] inputSlacks = new double[m];
System.arraycopy(solutionPoint, n + 1, inputSlacks, 0, m);
// 产出松弛变量
double[] outputSlacks = new double[s];
System.arraycopy(solutionPoint, n + 1 + m, outputSlacks, 0, s);
// lambda值
double[] lambdas = new double[n];
System.arraycopy(solutionPoint, 0, lambdas, 0, n);
return new DEAResult(theta, inputSlacks, outputSlacks, lambdas, dmuIndex, "BCC");
}
/**
* 计算所有DMU的BCC模型效率值
*/
public static List<DEAResult> calculateAllBCC(double[][] inputs, double[][] outputs) {
List<DEAResult> results = new ArrayList<>();
for (int i = 0; i < inputs.length; i++) {
results.add(calculateBCC(inputs, outputs, i));
}
return results;
}
// ===================== 非期望SBM模型 =====================
/**
* 计算指定DMU的非期望SBM模型效率值
* SBM模型考虑非期望产出如污染浪费等
*/
public static DEAResult calculateUndesirableSBM(
double[][] inputs,
double[][] desirableOutputs,
double[][] undesirableOutputs,
int dmuIndex) {
// 验证输入数据
validateSBMInput(inputs, desirableOutputs, undesirableOutputs, dmuIndex);
int n = inputs.length; // DMU数量
int m = inputs[0].length; // 投入指标数量
int s1 = desirableOutputs[0].length; // 期望产出指标数量
int s2 = undesirableOutputs[0].length; // 非期望产出指标数量
// SBM模型参数:lambda(n) + s^-(m) + s^+(s1) + s^b(s2)
int paramCount = n + m + s1 + s2;
// 初始值
double[] initialGuess = new double[paramCount];
for (int i = 0; i < paramCount; i++) {
initialGuess[i] = 0.1;
}
// 创建目标函数(显式实现MultivariateFunction接口,解决类型转换问题)
MultivariateFunction objectiveFunction = new SBMObjectiveFunction(
inputs, desirableOutputs, undesirableOutputs, dmuIndex, m, s1, s2, n);
// 创建优化器
MultivariateOptimizer optimizer = new SimplexOptimizer(1e-9, 1e-9);
NelderMeadSimplex simplex = new NelderMeadSimplex(paramCount, 0.1);
// 设置约束:参数非负
List<OptimizationData> optData = new ArrayList<>();
optData.add(new ObjectiveFunction(objectiveFunction));
optData.add(GoalType.MAXIMIZE); // SBM模型目标是最大化效率值
optData.add(simplex);
optData.add(new InitialGuess(initialGuess));
optData.add(new MaxIter(10000));
optData.add(new MaxEval(100000));
// 执行优化
PointValuePair solution = optimizer.optimize(optData.toArray(new OptimizationData[0]));
double[] solutionPoint = solution.getPoint();
double efficiency = solution.getValue();
// 解析结果
double[] lambda = new double[n];
double[] inputSlacks = new double[m];
double[] desirableOutputSlacks = new double[s1];
double[] undesirableOutputSlacks = new double[s2];
System.arraycopy(solutionPoint, 0, lambda, 0, n);
System.arraycopy(solutionPoint, n, inputSlacks, 0, m);
System.arraycopy(solutionPoint, n + m, desirableOutputSlacks, 0, s1);
System.arraycopy(solutionPoint, n + m + s1, undesirableOutputSlacks, 0, s2);
// 构建结果对象
DEAResult result = new DEAResult(efficiency, inputSlacks, desirableOutputSlacks,
lambda, dmuIndex, "UndesirableSBM");
result.setUndesirableOutputSlacks(undesirableOutputSlacks);
return result;
}
/**
* SBM模型的目标函数实现显式实现MultivariateFunction接口
*/
private static class SBMObjectiveFunction implements MultivariateFunction {
private final double[][] inputs;
private final double[][] desirableOutputs;
private final double[][] undesirableOutputs;
private final int dmuIndex;
private final int m;
private final int s1;
private final int s2;
private final int n;
public SBMObjectiveFunction(double[][] inputs, double[][] desirableOutputs,
double[][] undesirableOutputs, int dmuIndex,
int m, int s1, int s2, int n) {
this.inputs = inputs;
this.desirableOutputs = desirableOutputs;
this.undesirableOutputs = undesirableOutputs;
this.dmuIndex = dmuIndex;
this.m = m;
this.s1 = s1;
this.s2 = s2;
this.n = n;
}
@Override
public double value(double[] params) {
// 分离参数
double[] lambda = new double[n];
double[] sMinus = new double[m];
double[] sPlus = new double[s1];
double[] sBad = new double[s2];
System.arraycopy(params, 0, lambda, 0, n);
System.arraycopy(params, n, sMinus, 0, m);
System.arraycopy(params, n + m, sPlus, 0, s1);
System.arraycopy(params, n + m + s1, sBad, 0, s2);
// 计算分子:1 - (1/m)sum(s_i^-/x_ik)
double sumInputSlack = 0.0;
for (int i = 0; i < m; i++) {
sumInputSlack += sMinus[i] / inputs[dmuIndex][i];
}
double numerator = 1 - (sumInputSlack / m);
// 计算分母:1 + (1/(s1+s2))(sum(s_r^+/y_rk) + sum(s_b^b/z_bk))
double sumDesirableSlack = 0.0;
for (int r = 0; r < s1; r++) {
sumDesirableSlack += sPlus[r] / desirableOutputs[dmuIndex][r];
}
double sumUndesirableSlack = 0.0;
for (int b = 0; b < s2; b++) {
sumUndesirableSlack += sBad[b] / undesirableOutputs[dmuIndex][b];
}
double denominator = 1 + (sumDesirableSlack + sumUndesirableSlack) / (s1 + s2);
// 防止除零错误
if (denominator < 1e-10) {
return Double.NEGATIVE_INFINITY;
}
// 计算约束违反惩罚
double penalty = calculateConstraintPenalty(params, lambda, sMinus, sPlus, sBad);
// 返回目标函数值减去惩罚项
return (numerator / denominator) - penalty;
}
/**
* 计算约束违反的惩罚值
*/
private double calculateConstraintPenalty(double[] params, double[] lambda,
double[] sMinus, double[] sPlus, double[] sBad) {
double penalty = 0.0;
// 投入约束:sum(lambda_j * x_ij) + s_i^- = x_ik
for (int i = 0; i < m; i++) {
double sum = 0.0;
for (int j = 0; j < n; j++) {
sum += lambda[j] * inputs[j][i];
}
double constraint = sum + sMinus[i] - inputs[dmuIndex][i];
penalty += Math.abs(constraint) * 1e6; // 对违反约束进行惩罚
}
// 期望产出约束:sum(lambda_j * y_rj) - s_r^+ = y_rk
for (int r = 0; r < s1; r++) {
double sum = 0.0;
for (int j = 0; j < n; j++) {
sum += lambda[j] * desirableOutputs[j][r];
}
double constraint = sum - sPlus[r] - desirableOutputs[dmuIndex][r];
penalty += Math.abs(constraint) * 1e6;
}
// 非期望产出约束:sum(lambda_j * z_bj) + s_b^b = z_bk
for (int b = 0; b < s2; b++) {
double sum = 0.0;
for (int j = 0; j < n; j++) {
sum += lambda[j] * undesirableOutputs[j][b];
}
double constraint = sum + sBad[b] - undesirableOutputs[dmuIndex][b];
penalty += Math.abs(constraint) * 1e6;
}
// 变量非负约束惩罚
for (double param : params) {
if (param < 0) {
penalty += Math.abs(param) * 1e6;
}
}
// SBM模型的凸性约束 sum(lambda_j) = 1
double sumLambda = 0.0;
for (double l : lambda) {
sumLambda += l;
}
penalty += Math.abs(sumLambda - 1.0) * 1e6;
return penalty;
}
}
/**
* 计算所有DMU的非期望SBM模型效率值
*/
public static List<DEAResult> calculateAllUndesirableSBM(
double[][] inputs,
double[][] desirableOutputs,
double[][] undesirableOutputs) {
List<DEAResult> results = new ArrayList<>();
for (int i = 0; i < inputs.length; i++) {
results.add(calculateUndesirableSBM(inputs, desirableOutputs, undesirableOutputs, i));
}
return results;
}
// ===================== 验证方法 =====================
/**
* 验证SBM模型输入数据的合法性
*/
private static void validateSBMInput(double[][] inputs, double[][] desirableOutputs,
double[][] undesirableOutputs, int dmuIndex) {
validateInput(inputs, desirableOutputs, dmuIndex);
if (undesirableOutputs == null) {
throw new IllegalArgumentException("非期望产出矩阵不能为null");
}
if (undesirableOutputs.length != inputs.length) {
throw new IllegalArgumentException("非期望产出矩阵的行数必须与投入矩阵一致");
}
int s2 = undesirableOutputs[0].length;
for (double[] output : undesirableOutputs) {
if (output.length != s2) {
throw new IllegalArgumentException("所有DMU的非期望产出指标数量必须一致");
}
}
// 检查非期望产出是否为正值
for (int i = 0; i < undesirableOutputs.length; i++) {
for (int j = 0; j < s2; j++) {
if (undesirableOutputs[i][j] <= 0) {
throw new IllegalArgumentException("非期望产出必须为正值");
}
}
}
}
/**
* 验证基本DEA输入数据的合法性
*/
private static void validateInput(double[][] inputs, double[][] outputs, int dmuIndex) {
if (inputs == null || outputs == null) {
throw new IllegalArgumentException("投入和产出矩阵不能为null");
}
if (inputs.length != outputs.length) {
throw new IllegalArgumentException("投入和产出矩阵的行数(DMU数量)必须一致");
}
if (inputs.length == 0) {
throw new IllegalArgumentException("至少需要一个DMU");
}
int m = inputs[0].length;
for (double[] input : inputs) {
if (input.length != m) {
throw new IllegalArgumentException("所有DMU的投入指标数量必须一致");
}
for (double val : input) {
if (val <= 0) {
throw new IllegalArgumentException("投入必须为正值");
}
}
}
int s = outputs[0].length;
for (double[] output : outputs) {
if (output.length != s) {
throw new IllegalArgumentException("所有DMU的产出指标数量必须一致");
}
for (double val : output) {
if (val <= 0) {
throw new IllegalArgumentException("产出必须为正值");
}
}
}
if (dmuIndex < 0 || dmuIndex >= inputs.length) {
throw new IllegalArgumentException("DMU索引超出范围");
}
}
// ===================== 结果封装类 =====================
/**
* DEA计算结果封装类支持BCC和非期望SBM模型
*/
public static class DEAResult {
private final double efficiency; // 效率值
private final double[] inputSlacks; // 投入松弛变量
private final double[] outputSlacks; // 产出松弛变量(期望产出)
private double[] undesirableOutputSlacks; // 非期望产出松弛变量
private final double[] lambdas; // 各DMU的权重
private final int dmuIndex; // 对应的DMU索引
private final String modelType; // 模型类型
public DEAResult(double efficiency, double[] inputSlacks, double[] outputSlacks,
double[] lambdas, int dmuIndex, String modelType) {
this.efficiency = efficiency;
this.inputSlacks = inputSlacks;
this.outputSlacks = outputSlacks;
this.lambdas = lambdas;
this.dmuIndex = dmuIndex;
this.modelType = modelType;
}
/**
* 判断该DMU是否为DEA有效
*/
public boolean isEfficient() {
if (Math.abs(efficiency - 1.0) > 1e-9) {
return false;
}
// 检查所有松弛变量是否为0
for (double slack : inputSlacks) {
if (Math.abs(slack) > 1e-9) {
return false;
}
}
for (double slack : outputSlacks) {
if (Math.abs(slack) > 1e-9) {
return false;
}
}
// 检查非期望产出松弛变量
if (undesirableOutputSlacks != null) {
for (double slack : undesirableOutputSlacks) {
if (Math.abs(slack) > 1e-9) {
return false;
}
}
}
return true;
}
// setter和getter方法
public double getEfficiency() {
return efficiency;
}
public double[] getInputSlacks() {
return inputSlacks.clone();
}
public double[] getOutputSlacks() {
return outputSlacks.clone();
}
public double[] getUndesirableOutputSlacks() {
return undesirableOutputSlacks != null ? undesirableOutputSlacks.clone() : null;
}
public void setUndesirableOutputSlacks(double[] undesirableOutputSlacks) {
this.undesirableOutputSlacks = undesirableOutputSlacks;
}
public double[] getLambdas() {
return lambdas.clone();
}
public int getDmuIndex() {
return dmuIndex;
}
public String getModelType() {
return modelType;
}
@Override
public String toString() {
StringBuilder sb = new StringBuilder();
sb.append(modelType).append("模型 - DMU ").append(dmuIndex).append(" 评估结果:\n");
sb.append(" 效率值: ").append(String.format("%.4f", efficiency)).append("\n");
sb.append(" 是否DEA有效: ").append(isEfficient()).append("\n");
sb.append(" 投入松弛变量: ");
for (double slack : inputSlacks) {
sb.append(String.format("%.4f ", slack));
}
sb.append("\n 期望产出松弛变量: ");
for (double slack : outputSlacks) {
sb.append(String.format("%.4f ", slack));
}
if (undesirableOutputSlacks != null) {
sb.append("\n 非期望产出松弛变量: ");
for (double slack : undesirableOutputSlacks) {
sb.append(String.format("%.4f ", slack));
}
}
return sb.toString();
}
}
}
Loading…
Cancel
Save