1. 主程序框架
%% 基于遗传算法的债券投资组合优化
clear; clc; close all;
%% 债券数据生成（实际应用中替换为真实数据）
nBonds = 50; % 债券数目
bondData = generateBondData(nBonds);
%% 遗传算法参数设置
gaOptions = setupGAOptions();
%% 运行遗传算法优化
fprintf('开始债券投资组合优化...\n');
[optimalWeights, optimalFitness, exitFlag, output] = ...
runBondPortfolioGA(bondData, gaOptions);
%% 结果分析与可视化
analyzeResults(optimalWeights, optimalFitness, bondData, output);

2. 债券数据结构
function bondData = generateBondData(nBonds)
% 生成模拟债券数据
rng(42); % 设置随机种子保证可重复性
bondData = struct();
bondData.nBonds = nBonds;
% 债券基本属性
bondData.yield = 0.02 + 0.08 * rand(nBonds, 1); % 收益率 2%-10%
bondData.duration = 1 + 9 * rand(nBonds, 1); % 久期 1-10年
bondData.convexity = 0.5 + 4.5 * rand(nBonds, 1); % 凸性
bondData.creditRating = randi([1, 10], nBonds, 1); % 信用评级 (1最好,10最差)
bondData.sector = randi([1, 5], nBonds, 1); % 行业分类
bondData.liquidity = 0.1 + 0.9 * rand(nBonds, 1); % 流动性指标
% 风险参数
bondData.targetDuration = 5.0; % 目标久期
bondData.maxSectorWeight = 0.3; % 单一行业最大权重
bondData.budget = 1.0; % 总预算
end

3. 遗传算法配置
function gaOptions = setupGAOptions()
% 配置遗传算法参数
gaOptions = optimoptions('ga');
% 算法参数
gaOptions.PopulationSize = 100; % 种群规模
gaOptions.MaxGenerations = 200; % 最大代数
gaOptions.CrossoverFraction = 0.8; % 交叉比例
gaOptions.EliteCount = 2; % 精英个体数量
% 显示设置
gaOptions.Display = 'iter';
gaOptions.PlotFcn = {@gaplotbestf, @gaplotdistance};
% 停止条件
gaOptions.FunctionTolerance = 1e-6;
gaOptions.StallGenLimit = 50;
end

4. 适应度函数
function fitness = bondPortfolioFitness(weights, bondData)
% 债券投资组合适应度函数
% 确保权重为行向量
weights = weights(:)';
% 约束检查
[valid, penalty] = checkConstraints(weights, bondData);
if ~valid
fitness = -1e6 * penalty; % 严重惩罚违反约束
return;
end
% 计算收益目标
portfolioReturn = sum(weights .* bondData.yield');
% 计算风险惩罚项
riskPenalty = calculateRiskPenalty(weights, bondData);
% 适应度函数：收益 - 风险调整
fitness = portfolioReturn - riskPenalty;
end

function [valid, penalty] = checkConstraints(weights, bondData)
% 检查投资组合约束
penalty = 0;
% 预算约束
weightSum = sum(weights);

if abs(weightSum - bondData.budget) > 1e-6
penalty = penalty + 100 * abs(weightSum - bondData.budget);
end
% 非负限制
negativeWeights = sum(weights(weights < 0));
if negativeWeights < 0
penalty = penalty + 1000 * abs(negativeWeights);
end
% 单一债券权重限制 (最大20%)
maxWeightViolation = max(0, max(weights) - 0.2);
penalty = penalty + 500 * maxWeightViolation;
% 行业分散限制
sectorPenalty = checkSectorConstraints(weights, bondData);
penalty = penalty + sectorPenalty;
valid = (penalty == 0);
end

function sectorPenalty = checkSectorConstraints(weights, bondData)
% 检查行业分散限制
sectorPenalty = 0;
uniqueSectors = unique(bondData.sector);
for i = 1:length(uniqueSectors)
sectorMask = (bondData.sector == uniqueSectors(i));
sectorWeight = sum(weights(sectorMask));
% 检查是否超过最大行业权重
if sectorWeight > bondData.maxSectorWeight
sectorPenalty = sectorPenalty + 200 * (sectorWeight - bondData.maxSectorWeight);
end
end
end

function riskPenalty = calculateRiskPenalty(weights, bondData)
% 计算风险惩罚项
% 久期偏差惩罚
portfolioDuration = sum(weights .* bondData.duration');
durationDeviation = abs(portfolioDuration - bondData.targetDuration);
durationPenalty = 0.1 * durationDeviation;
% 信用风险惩罚
creditRisk = sum(weights .* (bondData.creditRating' / 10));
creditPenalty = 0.05 * creditRisk;
% 流动性惩罚
liquidityPenalty = 0.02 * sum(weights .* (1 - bondData.liquidity'));
% 集中度风险惩罚
concentrationPenalty = 0.1 * (sum(weights.^2));  % Herfindahl指数
riskPenalty = durationPenalty + creditPenalty + liquidityPenalty + concentrationPenalty;
end

5. 遗传算法执行
function [optimalWeights, optimalFitness, exitFlag, output] = runBondPortfolioGA(bondData, gaOptions)
% 运行遗传算法优化
nVars = bondData.nBonds;
% 约束定义
A = []; b = [];  % 线性不等式约束 A*x <= b
Aeq = ones(1, nVars);  % 权重和为1
beq = bondData.budget;
% 变量边界 [0, 0.2]
lb = zeros(1, nVars);
ub = 0.2 * ones(1, nVars);
% 非线性约束
nonlcon = @(weights) nonLinearConstraints(weights, bondData);
% 适应度函数
fitnessFcn = @(weights) -bondPortfolioFitness(weights, bondData);
% 运行遗传算法
[optimalWeights, optimalFitness, exitFlag, output] = ...

ga(fitnessFcn, nVars, A, b, Aeq, beq, lb, ub, nonlcon, gaOptions);
% 转换回正向适应度值
optimalFitness = -optimalFitness;
end

function [c, ceq] = nonLinearConstraints(weights, bondData)
% 非线性约束函数
c = [];  % 非线性不等式限制 c(x) <= 0
% 行业权重限制 (非线性等式限制)
uniqueSectors = unique(bondData.sector);
ceq = zeros(length(uniqueSectors), 1);
for i = 1:length(uniqueSectors)
    sectorMask = (bondData.sector == uniqueSectors(i));
    sectorWeight = sum(weights(sectorMask));
    % 允许轻微差异
    ceq(i) = max(0, sectorWeight - bondData.maxSectorWeight);
end
end

6. 结果分析
function analyzeResults(optimalWeights, optimalFitness, bondData, output)
% 分析并展示优化结果
fprintf('\n=== 优化结果分析 ===\n');
fprintf('最佳适应度值: %.4f\n', optimalFitness);
fprintf('投资组合收益率: %.4f%%\n', 100 * sum(optimalWeights .* bondData.yield'));
fprintf('投资组合久期: %.4f 年\n', sum(optimalWeights .* bondData.duration'));
% 投资组合统计
nSelected = sum(optimalWeights > 0.001);
fprintf('选定的债券数量: %d/%d\n', nSelected, bondData.nBonds);
% 权重分布
figure('Position', [100, 100, 1200, 800]);
subplot(2, 3, 1);
bar(optimalWeights(optimalWeights > 0.001));
title('债券权重分布');
xlabel('债券编号');
ylabel('权重');
grid on;
% 收益率 vs 久期散点图
subplot(2, 3, 2);
scatter(bondData.duration, bondData.yield, 50, optimalWeights, 'filled');
colorbar;
xlabel('久期 (年)');
ylabel('收益率');
title('债券特征与权重');
grid on;
% 行业分布
subplot(2, 3, 3);
uniqueSectors = unique(bondData.sector);
sectorWeights = zeros(size(uniqueSectors));
for i = 1:length(uniqueSectors)
    sectorMask = (bondData.sector == uniqueSectors(i));
    sectorWeights(i) = sum(optimalWeights(sectorMask));
end
pie(sectorWeights);
title('行业权重分布');
legend(arrayfun(@(x) sprintf('行业 %d', x), uniqueSectors, 'UniformOutput', false));
% 收敛曲线
subplot(2, 3, 4);
plot(output.fitness);
title('遗传算法收敛曲线');
xlabel('代数');
ylabel('最佳适应度');
grid on;
% 风险贡献分析
subplot(2, 3, 5);
riskContribution = optimalWeights .* bondData.creditRating';
riskContribution = riskContribution / sum(riskContribution);
bar(riskContribution(optimalWeights > 0.001));
title('信用风险贡献');
xlabel('债券编号');
ylabel('风险贡献比例');
grid on;
% 权重与流动性关系
subplot(2, 3, 6);

scatter(bondData.liquidity(optimalWeights > 0.001), ...
optimalWeights(optimalWeights > 0.001), 50, 'filled');
xlabel('流动性指标');
ylabel('权重');
title('权重 vs 流动性');
grid on;
% 详细统计信息
fprintf('\n投资组合详细统计:\n');
fprintf('信用风险加权均值: %.4f\n', sum(optimalWeights .* bondData.creditRating'));
fprintf('流动性加权均值: %.4f\n', sum(optimalWeights .* bondData.liquidity'));
fprintf('凸性加权均值: %.4f\n', sum(optimalWeights .* bondData.convexity'));
% 显示前10大持仓
[sortedWeights, sortedIdx] = sort(optimalWeights, 'descend');
fprintf('\n前10大持仓:\n');
fprintf('债券#\t权重\t收益率\t久期\t信用评级\n');
for i = 1:min(10, nSelected)
idx = sortedIdx(i);
fprintf('%d\t%.4f\t%.4f\t%.2f\t%d\n', ...
idx, sortedWeights(i), bondData.yield(idx), ...
bondData.duration(idx), bondData.creditRating(idx));
end
end
7. 高级功能扩展
%% 多目标优化版本
function [weights, scores] = multiObjectiveBondOptimization(bondData)
% 多目标优化：同时优化收益和风险
fitnessFcn = @(weights) [ -sum(weights .* bondData.yield'), ...  % 目标1: 最大化收益
calculateTotalRisk(weights, bondData) ]; % 目标2: 减少风险
nVars = bondData.nBonds;
A = []; b = [];
Aeq = ones(1, nVars);
beq = 1;
lb = zeros(1, nVars);
ub = 0.2 * ones(1, nVars);
[weights, scores] = gamultiobj(fitnessFcn, nVars, A, b, Aeq, beq, lb, ub);
end
function totalRisk = calculateTotalRisk(weights, bondData)
% 计算总风险度量
durationRisk = abs(sum(weights .* bondData.duration') - bondData.targetDuration);
creditRisk = sum(weights .* (bondData.creditRating' / 10));
concentrationRisk = sum(weights.^2);
totalRisk = 0.1 * durationRisk + 0.05 * creditRisk + 0.1 * concentrationRisk;
end