# Can not use gamma function in an external function

조회 수: 6 (최근 30일)
Geovane Gomes 2023년 9월 20일
댓글: Geovane Gomes 2023년 9월 20일
Dear,
I am trying to generate some random numbers with a specific distribution, but to calculate the parameters of this distribution I have to use the gamma function.
It works well when I use it in the same file, for example in live script, but from an external function does not work (see error below).
clear
mu = 0.15;
sigma = 0.03;
N = 1e3;
x0 = [2.001,1.0e3];
fun = @(x) sqrt(gamma(1-2/x)-(gamma(1-1/x)).^2)./gamma(1-1/x)-sigma/mu;
par2 = fzero(fun,x0);
par1 = mu/gamma(1-1/par2);
randomNumbers = gevrnd(1/par2,par1/par2,par1,N,1);
% checking
histogram(randomNumbers)
mu = mean(randomNumbers)
mu = 0.1482
sigma = std(randomNumbers)
sigma = 0.0275
Now trying to use it in an external function, I got the error 'Unrecognized function or variable 'gamma'.
Can someone explain why this happens?
##### 댓글 수: 2없음 표시없음 숨기기
Dyuman Joshi 2023년 9월 20일
"Now trying to use it in an external function, I got the error 'Unrecognized function or variable 'gamma'."
What is the external function? How are you calling it?
Geovane Gomes 2023년 9월 20일
This one:
obj.means = 0.15;
obj.standardDeviations = 0.03;
obj.distributions = "Frechet";
obj.nSimulations = 2e4;
randomNumbers = generaterandomnumbers(obj);
Error using fzero>localFirstFcnEvalError
FZERO cannot continue because user-supplied function_handle ==> @(x)sqrt(gamma(1-2/x)-(gamma(1-1/x)).^2)./gamma(1-1/x)-sigma/mu failed with the error below.

Unrecognized function or variable 'gamma'.

Error in fzero (line 231)
localFirstFcnEvalError(FunFcn,FunFcnIn,ME);

Error in solution>generaterandomnumbers (line 46)
par2 = fzero(fun,x0);
function randomNumbers = generaterandomnumbers(obj)
nVariables = length(obj.means);
% Pre-allocate memory for the random numbers
randomNumbers = zeros(obj.nSimulations, nVariables);
for i = 1:nVariables
switch lower(obj.distributions{i})
case 'normal'
randomNumbers(:,i) = normrnd(obj.means(i), ...
obj.standardDeviations(i), obj.nSimulations, 1);
case 'binomial'
randomNumbers(:,i) = binornd(obj.means(i), ...
obj.standardDeviations(i), obj.nSimulations, 1);
case 'lognormal'
mu = log((obj.means(i) ^ 2) / ...
sqrt(obj.standardDeviations(i) ^ 2 + obj.means(i) ^ 2));
sigma = sqrt(log(obj.standardDeviations(i) ^ 2 / ...
(obj.means(i) ^ 2) + 1));
randomNumbers(:,i) = lognrnd(mu, sigma, obj.nSimulations, 1);
case {'gumbel', 'extreme type 1'}
gamma = 0.577216;
sigmaHat = sqrt(6) * obj.standardDeviations(i) / pi;
muHat = obj.means(i) - gamma * sigmaHat;
pd = makedist('GeneralizedExtremeValue', 'k', 0, ...
'sigma', sigmaHat, 'mu', muHat);
randomNumbers(:,i) = random(pd, obj.nSimulations, 1);
case {'gumbel min', 'extreme type 1 min'}
gamma = 0.577216;
sigmaHat = sqrt(6) * obj.standardDeviations(i) / pi;
muHat = obj.means(i) + gamma / sigmaHat;
randomNumbers(:,i) = evrnd(muHat, sigmaHat, ...
obj.nSimulations, 1);
case 'frechet'
x0 = [2.001,1.0e3];
mu = obj.means(i);
sigma = obj.standardDeviations(i);
fun = @(x) sqrt(gamma(1-2/x)-(gamma(1-1/x)).^2)./gamma(1-1/x)-sigma/mu;
par2 = fzero(fun,x0);
par1 = mu/gamma(1-1/par2);
randomNumbers = gevrnd(1/par2,par1/par2,par1,N,1);
otherwise
error('Invalid distribution: %s', obj.distributions{i});
end
end
end

댓글을 달려면 로그인하십시오.

### 채택된 답변

Dyuman Joshi 2023년 9월 20일
편집: Dyuman Joshi 2023년 9월 20일
Using function names as variable names is not a good idea - in your case gamma.
I have changed the names to gamma0 and added a random value for N in the frechet case (as it was not defined), and the code runs without an errror.
obj.means = 0.15;
obj.standardDeviations = 0.03;
obj.distributions = "Frechet";
obj.nSimulations = 2e4;
randomNumbers = generaterandomnumbers(obj)
randomNumbers = 5×1
0.1178 0.1207 0.1649 0.1623 0.1244
function randomNumbers = generaterandomnumbers(obj)
nVariables = length(obj.means);
% Pre-allocate memory for the random numbers
randomNumbers = zeros(obj.nSimulations, nVariables);
for i = 1:nVariables
switch lower(obj.distributions{i})
case 'normal'
randomNumbers(:,i) = normrnd(obj.means(i), ...
obj.standardDeviations(i), obj.nSimulations, 1);
case 'binomial'
randomNumbers(:,i) = binornd(obj.means(i), ...
obj.standardDeviations(i), obj.nSimulations, 1);
case 'lognormal'
mu = log((obj.means(i) ^ 2) / ...
sqrt(obj.standardDeviations(i) ^ 2 + obj.means(i) ^ 2));
sigma = sqrt(log(obj.standardDeviations(i) ^ 2 / ...
(obj.means(i) ^ 2) + 1));
randomNumbers(:,i) = lognrnd(mu, sigma, obj.nSimulations, 1);
case {'gumbel', 'extreme type 1'}
%% Name changed
gamma0 = 0.577216;
sigmaHat = sqrt(6) * obj.standardDeviations(i) / pi;
muHat = obj.means(i) - gamma0 * sigmaHat;
pd = makedist('GeneralizedExtremeValue', 'k', 0, ...
'sigma', sigmaHat, 'mu', muHat);
randomNumbers(:,i) = random(pd, obj.nSimulations, 1);
case {'gumbel min', 'extreme type 1 min'}
%% Name changed
gamma0 = 0.577216;
sigmaHat = sqrt(6) * obj.standardDeviations(i) / pi;
muHat = obj.means(i) + gamma0 / sigmaHat;
randomNumbers(:,i) = evrnd(muHat, sigmaHat, ...
obj.nSimulations, 1);
case 'frechet'
x0 = [2.001,1.0e3];
mu = obj.means(i);
sigma = obj.standardDeviations(i);
fun = @(x) sqrt(gamma(1-2/x)-(gamma(1-1/x)).^2)./gamma(1-1/x)-sigma/mu;
par2 = fzero(fun,x0);
par1 = mu/gamma(1-1/par2);
N=5;
randomNumbers = gevrnd(1/par2,par1/par2,par1,N,1);
otherwise
error('Invalid distribution: %s', obj.distributions{i});
end
end
end
##### 댓글 수: 2없음 표시없음 숨기기
Dyuman Joshi 2023년 9월 20일
Though, variables defined within a case are not accessible to another case and MATLAB executes only one case of any switch statement.
So for the input you have provided, only the 2nd last case is executed and gamma is not over-writed. But the error still occurs. Hmmm.
Geovane Gomes 2023년 9월 20일
Thank you very much for your help. It works.
"N" is actually "obj.nSimulations". It was wrong. Thanks again.

댓글을 달려면 로그인하십시오.

### 카테고리

Help CenterFile Exchange에서 Gamma Distribution에 대해 자세히 알아보기

R2023b

### Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by