How to graph R0 against two parameters on a 3D plot with planes for certain R0 values?

I have been trying to replicate the following graph from the paper, https://pubmed.ncbi.nlm.nih.gov/32834653/. The authors have plotted R0 against the two parameters gamma_a and beta_s on a 3D plot along with two planes for the values of R0=1 and R0=2. The parameters values and the formula for R0 are given below.

beta_s = 0.274;
alpha_a = 0.4775;
alpha_u = 0.695;
mu = 0.062;
q_i = 0.078;
gamma_a = 0.29;
eta_i = 0.009;
eta_u = 0.05;

R0 = (beta_s*alpha_a)/(gamma_a+mu) + (beta_s*alpha_u*gamma_a*(1-q_i))/((gamma_a+mu)*(eta_u+mu))

I would highly appreciate if someone could help me draw this plot. Thank you so much!

Answer 1

[Bora Eryilmaz]

[X,Y] = meshgrid(0:0.1:0.4, 0:0.01:0.4);
R0 = 1;
Z = R0 * ones(size(X));
surf(X,Y,Z)
colormap('turbo')
shading interp
xlabel('\beta_s')
ylabel('\gamma_a')
zlabel('R_0')
hold on
    
R0 = 2;
Z = R0 * ones(size(X));
surf(X,Y,Z)
colormap('turbo')
shading interp

%beta_s = 0.274;
beta_s = X;
alpha_a = 0.4775;
alpha_u = 0.695;
mu = 0.062;
q_i = 0.078;
% gamma_a = 0.29;
gamma_a = Y;
eta_i = 0.009;
eta_u = 0.05;

R0 = (beta_s.*alpha_a)./(gamma_a+mu) + (beta_s.*alpha_u.*gamma_a*(1-q_i))./((gamma_a+mu)*(eta_u+mu));
surf(X,Y,R0)
colormap('turbo')
shading interp

[Bas123]

It worked great. I do sincerely appreciate your support. Thanks a lot!