Skip to main content
added 1792 characters in body
Source Link
Sextus Empiricus
  • 86.5k
  • 6
  • 115
  • 304

A gamma distribution with given median and IQR

f = function(par, targetmedian, targetIQR) {
  k = par[1]
  theta = par[2]

  median = theta*k*(1-1/9/k)^3
  ### alternatively compute the median like below
  #median = qgamma(0.5,shape = k,scale = theta)

  IQR = qgamma(0.75,shape = k,scale = theta) - qgamma(0.25,shape = k,scale = theta)

  error = (median-targetmedian)^2+(IQR-targetIQR)^2
  return(error)
}

### optim function below gives result k = 2.905083 theta = 1.160847
p2 = optim(par  = c(1,1), f, targetmedian = 3, targetIQR = 2.5)$par

### check median and IQR
qgamma(0.5, shape = p2[1], scale = p2[2]) ### 2.994302
qgamma(0.75,shape = p2[1], scale = p2[2])-qgamma(0.25,shape = p1[1],scale = p1[2]) # 2.494646

A gamma distribution with given quartiles

In this case you do not have two, but three values, whereas the gamma distribution only has two parameters. It might be possible that an exact fit is not possible. Now the cost function for the optimization becomes more important.

  • You could use the Kolmogorov statistic (the distance between the cdf and your quartiles that define a empirical distribution)

  • You can consider data bins $x\leq 2$, $2<x\leq 3$, $3 < x \leq 5$ and $5<x$ with each 25% probability and optimize the cross-entropy.

example code:

f = function(par, target_quartiles, fit_method = "Kolmogorov") {
  k = par[1]
  theta = par[2]

  if (fit_method %in% c("Kolmogorov", "cross-entropy") == FALSE) { 
    stop("use method Kolmogorov or cross-entropy")
  }

  p_quartiles = pgamma(target_quartiles, shape = k, scale = theta)
  
  if (fit_method == "Kolmogorov") {
     statistic = max(abs(p_quartiles - c(0.25,0.5,0.75)))
  }
  if (fit_method == "cross-entropy") {
     p1 = diff(c(0,p_quartiles,1))
     p2 = c(0.25,0.25,0.25,0.25)
     statistic = -sum(p2*log(p1))
  }     
  return(statistic)
}   

### optim function below gives result k = 2.399316 theta = 1.502534
par1 = optim(par  = c(1,1), f, target_quartiles = c(2,3,5))$par

### 0.2724194 0.4775805 0.7724194
pgamma(c(2,3,5), shape = par1[1], scale = par1[2])

### optim function below gives result k = 2.408474 theta = 1.530167
par2 = optim(par  = c(1,1), f, target_quartiles = c(2,3,5), fit_method = "cross-entropy")$par

### 0.2625420 0.4645755 0.7610788
pgamma(c(2,3,5), shape = par2[1], scale = par2[2])
f = function(par, targetmedian, targetIQR) {
  k = par[1]
  theta = par[2]

  median = theta*k*(1-1/9/k)^3
  ### alternatively compute the median like below
  #median = qgamma(0.5,shape = k,scale = theta)

  IQR = qgamma(0.75,shape = k,scale = theta) - qgamma(0.25,shape = k,scale = theta)

  error = (median-targetmedian)^2+(IQR-targetIQR)^2
  return(error)
}

### optim function below gives result k = 2.905083 theta = 1.160847
p2 = optim(par  = c(1,1), f, targetmedian = 3, targetIQR = 2.5)$par

### check median and IQR
qgamma(0.5, shape = p2[1], scale = p2[2]) ### 2.994302
qgamma(0.75,shape = p2[1], scale = p2[2])-qgamma(0.25,shape = p1[1],scale = p1[2]) # 2.494646

A gamma distribution with given median and IQR

f = function(par, targetmedian, targetIQR) {
  k = par[1]
  theta = par[2]

  median = theta*k*(1-1/9/k)^3
  ### alternatively compute the median like below
  #median = qgamma(0.5,shape = k,scale = theta)

  IQR = qgamma(0.75,shape = k,scale = theta) - qgamma(0.25,shape = k,scale = theta)

  error = (median-targetmedian)^2+(IQR-targetIQR)^2
  return(error)
}

### optim function below gives result k = 2.905083 theta = 1.160847
p2 = optim(par  = c(1,1), f, targetmedian = 3, targetIQR = 2.5)$par

### check median and IQR
qgamma(0.5, shape = p2[1], scale = p2[2]) ### 2.994302
qgamma(0.75,shape = p2[1], scale = p2[2])-qgamma(0.25,shape = p1[1],scale = p1[2]) # 2.494646

A gamma distribution with given quartiles

In this case you do not have two, but three values, whereas the gamma distribution only has two parameters. It might be possible that an exact fit is not possible. Now the cost function for the optimization becomes more important.

  • You could use the Kolmogorov statistic (the distance between the cdf and your quartiles that define a empirical distribution)

  • You can consider data bins $x\leq 2$, $2<x\leq 3$, $3 < x \leq 5$ and $5<x$ with each 25% probability and optimize the cross-entropy.

example code:

f = function(par, target_quartiles, fit_method = "Kolmogorov") {
  k = par[1]
  theta = par[2]

  if (fit_method %in% c("Kolmogorov", "cross-entropy") == FALSE) { 
    stop("use method Kolmogorov or cross-entropy")
  }

  p_quartiles = pgamma(target_quartiles, shape = k, scale = theta)
  
  if (fit_method == "Kolmogorov") {
     statistic = max(abs(p_quartiles - c(0.25,0.5,0.75)))
  }
  if (fit_method == "cross-entropy") {
     p1 = diff(c(0,p_quartiles,1))
     p2 = c(0.25,0.25,0.25,0.25)
     statistic = -sum(p2*log(p1))
  }     
  return(statistic)
}   

### optim function below gives result k = 2.399316 theta = 1.502534
par1 = optim(par  = c(1,1), f, target_quartiles = c(2,3,5))$par

### 0.2724194 0.4775805 0.7724194
pgamma(c(2,3,5), shape = par1[1], scale = par1[2])

### optim function below gives result k = 2.408474 theta = 1.530167
par2 = optim(par  = c(1,1), f, target_quartiles = c(2,3,5), fit_method = "cross-entropy")$par

### 0.2625420 0.4645755 0.7610788
pgamma(c(2,3,5), shape = par2[1], scale = par2[2])
edited body
Source Link
Sextus Empiricus
  • 86.5k
  • 6
  • 115
  • 304

It is not easy because there is no closed form expression for the median and quartiles. For the median you can find several approximate formula's. For the quartiles you could use R's approximation of the quantile function, qgamma, and then use an iterative algorithm that finds $k$ and $\theta$ that suit your target.

One way to make such algorithm can be by optimizing a cost function

f = function(par, targetmedian, targetIQR) {
  k = par[1]
  theta = par[2]

  median = theta*k*(1-1/9/k)^3
  ### alternatively compute the median like below
  #median = qgamma(0.5,shape = k,scale = theta)

  IQR = qgamma(0.75,shape = k,scale = theta) - qgamma(0.25,shape = k,scale = theta)

  error = (median-targetmedian)^2+(IQR-targetIQR)^2
  return(error)
}

### optim function below gives result k = 2.905083 theta = 1.160847
p2 = optim(par  = c(1,1), f, targetmedian = 53, targetIQR = 2.5)$par

### check median and IQR
qgamma(0.5, shape = p2[1], scale = p2[2]) ### 2.994302
qgamma(0.75,shape = p2[1], scale = p2[2])-qgamma(0.25,shape = p1[1],scale = p1[2]) # 2.494646

It is not easy because there is no closed form expression for the median and quartiles. For the median you can find several approximate formula's. For the quartiles you could use R's approximation of the quantile function, qgamma, and then use an iterative algorithm that finds $k$ and $\theta$ that suit your target.

One way to make such algorithm can be by optimizing a cost function

f = function(par, targetmedian, targetIQR) {
  k = par[1]
  theta = par[2]

  median = theta*k*(1-1/9/k)^3
  ### alternatively compute the median like below
  #median = qgamma(0.5,shape = k,scale = theta)

  IQR = qgamma(0.75,shape = k,scale = theta) - qgamma(0.25,shape = k,scale = theta)

  error = (median-targetmedian)^2+(IQR-targetIQR)^2
  return(error)
}

### optim function below gives result k = 2.905083 theta = 1.160847
p2 = optim(par  = c(1,1), f, targetmedian = 5, targetIQR = 2.5)$par

### check median and IQR
qgamma(0.5, shape = p2[1], scale = p2[2]) ### 2.994302
qgamma(0.75,shape = p2[1], scale = p2[2])-qgamma(0.25,shape = p1[1],scale = p1[2]) # 2.494646

It is not easy because there is no closed form expression for the median and quartiles. For the median you can find several approximate formula's. For the quartiles you could use R's approximation of the quantile function, qgamma, and then use an iterative algorithm that finds $k$ and $\theta$ that suit your target.

One way to make such algorithm can be by optimizing a cost function

f = function(par, targetmedian, targetIQR) {
  k = par[1]
  theta = par[2]

  median = theta*k*(1-1/9/k)^3
  ### alternatively compute the median like below
  #median = qgamma(0.5,shape = k,scale = theta)

  IQR = qgamma(0.75,shape = k,scale = theta) - qgamma(0.25,shape = k,scale = theta)

  error = (median-targetmedian)^2+(IQR-targetIQR)^2
  return(error)
}

### optim function below gives result k = 2.905083 theta = 1.160847
p2 = optim(par  = c(1,1), f, targetmedian = 3, targetIQR = 2.5)$par

### check median and IQR
qgamma(0.5, shape = p2[1], scale = p2[2]) ### 2.994302
qgamma(0.75,shape = p2[1], scale = p2[2])-qgamma(0.25,shape = p1[1],scale = p1[2]) # 2.494646
added 77 characters in body
Source Link
Sextus Empiricus
  • 86.5k
  • 6
  • 115
  • 304

It is not easy because there is no closed form expression for the median and quartiles. For the median you can find several approximate formula's. For the quartiles you could use R's approximation of the quantile function, qgamma, and then use an iterative algorithm that finds $k$ and $\theta$ that suit your target.

One way to make such algorithm can be by optimizing a cost function

f = function(par, targetmedian, targetIQR) {
  k = par[1]
  theta = par[2]

  median = theta*k*(1-1/9/k)^3
  ### alternatively compute the median like below
  #median = qgamma(0.5,shape = k,scale = theta)

  IQR = qgamma(0.75,shape = k,scale = theta) - qgamma(0.25,shape = k,scale = theta)

  error = (median-3targetmedian)^2+(IQR-2.5targetIQR)^2
  return(error)
}

### optim function below gives result k = 2.905083 theta = 1.160847
p2 = optim(par  = c(1,1), f, targetmedian = 5, targetIQR = 2.5)$par

### check median and IQR
qgamma(0.5, shape = p2[1], scale = p2[2]) ### 2.994302
qgamma(0.75,shape = p2[1], scale = p2[2])-qgamma(0.25,shape = p1[1],scale = p1[2]) # 2.494646

It is not easy because there is no closed form expression for the median and quartiles. For the median you can find several approximate formula's. For the quartiles you could use R's approximation of the quantile function, qgamma, and then use an iterative algorithm that finds $k$ and $\theta$ that suit your target.

One way to make such algorithm can be by optimizing a cost function

f = function(par) {
  k = par[1]
  theta = par[2]

  median = theta*k*(1-1/9/k)^3
  ### alternatively compute the median like below
  #median = qgamma(0.5,shape = k,scale = theta)

  IQR = qgamma(0.75,shape = k,scale = theta) - qgamma(0.25,shape = k,scale = theta)

  error = (median-3)^2+(IQR-2.5)^2
  return(error)
}

### result k = 2.905083 theta = 1.160847
p2 = optim(par  = c(1,1), f)$par
qgamma(0.5, shape = p2[1], scale = p2[2]) ### 2.994302
qgamma(0.75,shape = p2[1], scale = p2[2])-qgamma(0.25,shape = p1[1],scale = p1[2]) # 2.494646

It is not easy because there is no closed form expression for the median and quartiles. For the median you can find several approximate formula's. For the quartiles you could use R's approximation of the quantile function, qgamma, and then use an iterative algorithm that finds $k$ and $\theta$ that suit your target.

One way to make such algorithm can be by optimizing a cost function

f = function(par, targetmedian, targetIQR) {
  k = par[1]
  theta = par[2]

  median = theta*k*(1-1/9/k)^3
  ### alternatively compute the median like below
  #median = qgamma(0.5,shape = k,scale = theta)

  IQR = qgamma(0.75,shape = k,scale = theta) - qgamma(0.25,shape = k,scale = theta)

  error = (median-targetmedian)^2+(IQR-targetIQR)^2
  return(error)
}

### optim function below gives result k = 2.905083 theta = 1.160847
p2 = optim(par  = c(1,1), f, targetmedian = 5, targetIQR = 2.5)$par

### check median and IQR
qgamma(0.5, shape = p2[1], scale = p2[2]) ### 2.994302
qgamma(0.75,shape = p2[1], scale = p2[2])-qgamma(0.25,shape = p1[1],scale = p1[2]) # 2.494646
Source Link
Sextus Empiricus
  • 86.5k
  • 6
  • 115
  • 304
Loading