If you have a Hyperexponential-2 (H2), then with probability $p$ you sample from $F_1$ and with $1-p$ sample from $F_2$. Obviously $p$ must be on the interval $[0,1]$, $F_1$ ~ Exponential($\lambda_1$), and $F_2$ ~ Exponential($\lambda_2$).
You can mix a larger number of exponentials by extending this idea.
If you have a target mean and variance, select $p,\lambda_1,\lambda_2$ so that the resulting hyperexponential will have the target mean and variance.
To choose the parameters, solve the equations for the mean and variance. If $X$ ~ $H2(p,\lambda_1,\lambda_2)$ then $\text{E}[X] = p/\lambda_1 + (1-p)/\lambda_2$.
$\text{Var}[X] = \text{E}[X^2]-\text{E}[X]^2 = 2p/\lambda_1 ^2 + 2(1-p)/\lambda_2^2 - \text{E}[X]^2$
Set these equal to your targets and solve for $p,\lambda_1,\lambda_2$; it will be underdetermined so a number of hyperexponentials will have your target mean and variance. Just pick one solution.
MATLAB Code to generate from $H2(p,\lambda_1,\lambda_2)$:
function [ X ] = h2rnd( p,Rates,N )
%H2RAND Samples from specified Hyperexponential distrbution
% [ X ] = h2rnd(p,Rates,N)
% INPUTS
% p: mixing probability
% Rates: rates for the two exponential distributions (2 x 1 vector)
% N: number of samples to generate
% OUTPUTS
% X ~ H2(p,Rates(1),Rates(2))
% %%%%%%%%%%%%% BEGIN ERROR CHECKING %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
if p < 0 | p > 1, error('p must be on interval [0,1]'), end
if min(Rates)<=0, error('Rates must be positive'), end
% %%%%%%%%%%%%% END ERROR CHECKING %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
X = zeros(N,1);
P = rand(N,1);
X(P<=p)= -log(1-rand(sum(P<=p),1))/Rates(1);
X(P>p)= -log(1-rand(sum(P>p),1))/Rates(2);
end
Update:
Adding another exponential to have $Y\sim H3(p_1,p_2,\lambda_1,\lambda_2,\lambda_3)$ allows one to match more moments but increases the complexity as well.
Reference:
Wiki for Hyperexponential Distribution
rexp(n,lambda[sample(p,n,TRUE,p)])
generates n values from the hyperexponential if you've already specifiedp
,lambda
andn
. ] $\endgroup$