%% 分析非诚实收益N对演化的影响
function analyze_dishonest_gain_parameter()
    
    fprintf('\n开始非诚实收益N的敏感性分析...\n');
    
    % 基础参数设置
    base_params = struct();
    base_params.k1 = 1;
    base_params.C2 = 50;
    base_params.D = 60;
    base_params.C3 = 50;
    base_params.p1 = 200;
    base_params.C = 40;
    base_params.m = 100;
    base_params.a = 0.5;
    base_params.k2 = 1;
    
    % 两组N值设置
    N_groups = {
        [50, 70, 100],     % 第一组：较小的N值
        [100, 150, 200]    % 第二组：较大的N值
    };
    
    % 对应的n值设置
    n_values = [50, 50];  % 两组使用不同的n值
    
    % 创建图形窗口
    figure('Position', [100, 100, 1200, 600]);
    
    for group_idx = 1:2
        subplot(1, 2, group_idx);
        hold on;
        
        N_values = N_groups{group_idx};
        base_params.n = n_values(group_idx);
        
        % 颜色设置
        colors = {[0.2, 0.4, 0.8], [0.4, 0.6, 0.2], [0.8, 0.4, 0.2]};
        line_styles = {'-', '--', ':'};
        
        % 时间跨度
        tspan = [0, 50];
        
        for N_idx = 1:3
            N = N_values(N_idx);
            
            % 运行多次仿真以展示不同初始条件的影响
            num_sims = 10;
            
            for sim = 1:num_sims
                % 随机初始条件
                initial_state = [rand(); rand(); rand()];
                
                % 求解演化方程
                [t, trajectory] = ode45(@(t, state) evolutionary_dynamics_N(t, state, base_params, N), ...
                                        tspan, initial_state);
                
                % 绘制轨迹（使用透明度）
                if sim == 1
                    h = plot3(trajectory(:, 1), trajectory(:, 2), trajectory(:, 3), ...
                             'Color', colors{N_idx}, 'LineWidth', 2, ...
                             'LineStyle', line_styles{N_idx});
                else
                    plot3(trajectory(:, 1), trajectory(:, 2), trajectory(:, 3), ...
                         'Color', [colors{N_idx}, 0.3], 'LineWidth', 1, ...
                         'LineStyle', line_styles{N_idx});
                end
            end
        end
        
        % 设置图形属性
        xlabel('x (监管部门)');
        ylabel('y (共享方A)');
        zlabel('z (共享方B)');
        title(sprintf('第%d组: N的影响', group_idx));
        grid on;
        view(45, 30);
        xlim([0, 1]);
        ylim([0, 1]);
        zlim([0, 1]);
        
        % 添加图例
        legend(arrayfun(@(i) sprintf('N=%d', N_values(i)), 1:3, 'UniformOutput', false), ...
               'Location', 'best');
        
        hold off;
    end
    
    sgtitle('非诚实收益N的敏感性分析（图9）', 'FontSize', 14, 'FontWeight', 'bold');
    
    % 生成2D投影图
    generate_N_2D_analysis(base_params, N_groups, n_values);
end

function generate_N_2D_analysis(base_params, N_groups, n_values)
    figure('Position', [100, 100, 1200, 600]);
    
    for group_idx = 1:2
        subplot(1, 2, group_idx);
        hold on;
        
        N_values = N_groups{group_idx};
        base_params.n = n_values(group_idx);
        
        % 创建时间向量
        tspan = [0, 50];
        t_plot = linspace(0, 50, 100);
        
        % 颜色和标记设置
        colors = {'r', 'g', 'b'};
        markers = {'o', 's', 'd'};
        
        for N_idx = 1:3
            N = N_values(N_idx);
            
            % 固定初始条件以便比较
            initial_state = [0.5; 0.5; 0.5];
            
            % 求解演化方程
            [t, trajectory] = ode45(@(t, state) evolutionary_dynamics_N(t, state, base_params, N), ...
                                    tspan, initial_state);
            
            % 绘制各个策略随时间的演化
            plot(t, trajectory(:, 1), [colors{N_idx}, '-'], 'LineWidth', 2);
            plot(t, trajectory(:, 2), [colors{N_idx}, '--'], 'LineWidth', 2);
            plot(t, trajectory(:, 3), [colors{N_idx}, ':'], 'LineWidth', 2);
            
            % 标记终点
            scatter(t(end), trajectory(end, 1), 100, colors{N_idx}, markers{N_idx}, 'filled');
            scatter(t(end), trajectory(end, 2), 100, colors{N_idx}, markers{N_idx}, 'filled');
            scatter(t(end), trajectory(end, 3), 100, colors{N_idx}, markers{N_idx}, 'filled');
        end
        
        % 设置图形属性
        xlabel('时间 t');
        ylabel('策略概率');
        title(sprintf('第%d组: 策略演化时间序列', group_idx));
        grid on;
        xlim([0, 50]);
        ylim([0, 1]);
        
        % 添加图例
        legend_entries = {};
        for i = 1:3
            legend_entries{3*(i-1)+1} = sprintf('N=%d (x)', N_values(i));
            legend_entries{3*(i-1)+2} = sprintf('N=%d (y)', N_values(i));
            legend_entries{3*(i-1)+3} = sprintf('N=%d (z)', N_values(i));
        end
        legend(legend_entries, 'Location', 'best', 'NumColumns', 3);
        
        hold off;
    end
    
    sgtitle('非诚实收益N对策略演化的影响', 'FontSize', 14, 'FontWeight', 'bold');
end

function dstate_dt = evolutionary_dynamics_N(~, state, params, N)
    x = state(1);
    y = state(2);
    z = state(3);
    
    % 提取参数
    k1 = params.k1;
    k2 = params.k2;
    C2 = params.C2;
    C3 = params.C3;
    D = params.D;
    p1 = params.p1;
    C = params.C;
    m = params.m;
    a = params.a;
    n = params.n;
    
    % 复制动态方程（包含N的影响）
    dx_dt = x * (1 - x) * ((C2 * D * y * z) + D - D * y * z - C3 + C2 * z * y);
    dy_dt = y * (1 - y) * (k1 * p1 * z + C * z + m * x + a * N - N - C);
    dz_dt = z * (1 - z) * (a * N + C * y + m * x + k2 * p1 * y - N - C);
    
    dstate_dt = [dx_dt; dy_dt; dz_dt];
end