请求协助绘制不同β值下SIR模型的感染人群时间曲线

I see exactly what you need—you want to plot Time vs Infected population curves for different β values (like 0.1, 0.2, 0.4, etc.) on the same graph, instead of plotting all 5 population compartments at once. Let's adjust your code to make that happen.

Key Issues in Your Original Code #

• It only runs for a single hard-coded β value (beta = 0.9)
• It plots all 5 population variables, but you only care about the Infected group (Y(:,2))
• No mechanism to overlay multiple curves for different βs

Modified Working Code #

function EulerSIRTenga()
    clear all; clc;
    rho = 0.7; Lambda = 0.4; 
    theta1 = 0.1; theta2 = 0.2; theta3 = 0.1;
    alpha = 0.33; omega = 0.4; psi = 0.3;
    
    % Define the beta values you want to compare
    beta_values = [0.1, 0.2, 0.4, 0.9];
    % Define distinct line styles/colors for each beta (for clarity)
    line_styles = {'r-', 'b--', 'g-.', 'k:'};
    
    % Initialize figure and hold to overlay all curves
    figure('Position', [100 100 800 500]);
    hold on;
    grid on;
    box on;
    
    % Loop through each beta value
    for i = 1:length(beta_values)
        beta = beta_values(i);
        % Calculate R0 for current beta (optional, prints to console)
        R0 = rho * beta/(theta1 + Lambda);
        fprintf('For β = %.1f, R0 = %.4f\n', beta, R0);
        
        % Solve ODE with current beta
        options = odeset('RelTol',1e-4,'AbsTol',[1e-4 1e-4 1e-4 1e-4 1e-4]);
        [T,Y] = ode45(@(t,y) SIRmodel(t,y,beta), [0 365], [0.9,0.7,0.5,0.15,0.35], options);
        
        % Plot only Time vs Infected (Y(:,2)) with unique style
        plot(T, Y(:,2), line_styles{i}, 'Linewidth', 2);
    end
    
    % Add labels, title, and legend
    xlabel('Time (days)','FontSize',14);
    ylabel('Proportion of Infected Population','FontSize',14);
    title('Infected Population vs Time for Different β Values','FontSize',14);
    % Create legend using beta values
    legend_strs = arrayfun(@(b) sprintf('β = %.1f', b), beta_values, 'UniformOutput', false);
    legend(legend_strs, 'FontSize', 12);
    
    hold off;
end

function dy=SIRmodel(t,y,beta)
    dy=zeros(5,1);
    rho = 0.7; Lambda = 0.4; 
    theta1 = 0.1; theta2 = 0.2; theta3 = 0.1;
    alpha = 0.33; omega = 0.4; psi = 0.3;
    
    dy(1)= rho*Lambda - (beta*y(2)+Lambda-alpha*(y(3)+y(4)+y(5)))*y(1);
    dy(2)=(1-rho)*Lambda+(beta*y(1)-theta1-Lambda+alpha*(y(3)+y(4)+y(5)))*y(2);
    dy(3)= theta1*y(2)+ theta2*y(4)+theta3*y(5)-(alpha+omega+Lambda-alpha*(y(3)+y(4)+y(5)))*y(3);
    dy(4)= psi*omega*y(3)-(alpha+theta2+Lambda-alpha*(y(3)+y(4)+y(5)))*y(4);
    dy(5)=(1-psi)*omega*y(3)-(alpha-theta3+Lambda-alpha*(y(3)+y(4)+y(5)))*y(5);
end

What Changed? #

1. β Value Array: We defined beta_values with all the values you want to test—you can easily add/remove values here to expand your comparison.
2. Parameter Passing: Modified SIRmodel to accept beta as an input parameter, and used @(t,y) SIRmodel(t,y,beta) to pass the current β to the ODE solver.
3. Overlay Plotting: Used hold on to keep the figure open while plotting each β's curve, with unique line styles/colors to make each curve easy to distinguish.
4. Focused Output: Only plots the Infected population (Y(:,2)) instead of all 5 compartments, matching your requirement.
5. Clear Legend: Automatically generates a legend labeled with each β value for quick, intuitive comparison.

How to Use #

Just run the EulerSIRTenga function, and it will generate a single graph with all your β curves, showing how infection levels change over time for different transmission rates.

内容的提问来源于stack exchange，提问作者Peter Tengaa