Bayesian Inference of a Finite Population Mean Under Length-Biased   Sampling

Zhiqing Xu; Balgobin Nandram; Binod Manandhar

arXiv:1902.04242·stat.AP·February 13, 2019

Bayesian Inference of a Finite Population Mean Under Length-Biased Sampling

Zhiqing Xu, Balgobin Nandram, Binod Manandhar

PDF

Open Access

TL;DR

This paper introduces a Bayesian approach for estimating the average shrub width in a finite population with unknown size, using length-biased sampling data, and demonstrates its robustness and improved performance over non-length-biased models.

Contribution

The paper develops a Bayesian method that accounts for length-biased sampling and unknown population size, employing a generalized gamma distribution for robustness.

Findings

01

Model with length bias performs better than without.

02

The method effectively estimates total shrub area.

03

Robust Bayesian inference via a random sampler enhances computation.

Abstract

We present a robust Bayesian method to analyze forestry data when samples are selected with probability proportional to length from a finite population of unknown size. Specifically, we use Bayesian predictive inference to estimate the finite population mean of shrub widths in a limestone quarry dominated by re-growth of mountain mahogany. The data on shrub widths are collected using transect sampling and it is assumed that the probability that a shrub is selected is proportional to its width; this is length-biased sampling. In this type of sampling, the population size is also unknown and this creates an additional challenge. The quantity of interest is average finite population shrub width and the total shrub area of the quarry can be estimated. Our method is assisted by using the three-parameter generalized gamma distribution, thereby robustifying our procedure against a possible…

Tables6

Table 1. Table 1: Widths (meters) of shrubs in Replication 1

Transect	$X_{i}$ =width
I	1.53, 0.87, 0.79, 0.78, 1.85, 1.45
	0.48, 0.52, 0.22, 0.38, 0.59, 0.20
	0.42, 1.02, 0.97, 0.56, 0.62, 0.42
II	1.15, 0.87, 0.57, 0.97, 0.57, 1.97
	0.58, 2.54, 1.85, 0.35, 1.24, 1.80
	0.78, 0.98, 1.30, 1.55, 1.69, 2.12
	1.27, 0.75, 1.01, 1.82
III	0.71, 1.50, 1.82, 1.86, 1.61, 1.21

Table 2. Table 2: Widths (meters) of shrubs in Replication 2

Transect	$X_{i}$ =width
I	0.67, 0.31, 0.83, 1.95, 1.36, 1.45
	0.72, 1.15, 0.98, 1.29, 0.88, 0.25
	0.63, 1.12, 0.34, 0.21, 1.36, 0.95
	1.04, 0.48, 1.05, 0.88, 0.16, 1.08
	0.95, 0.25, 0.30, 1.40, 0.58, 0.73
	1.30, 0.57
II	0.96, 2.08, 0.68, 1.39, 0.50, 0.72
	0.19, 1.91, 0.88, 0.48, 0.12

Table 3. Table 3: Estimated parameters and population mean in the model with corrected selection bias. 95% Bootstrap Confidence Intervals of x ¯ ¯ 𝑥 \bar{x} is (.86, 1.11).

$α$	$β_{1}$	$β_{2}$	$β_{3}$	$γ$	$\bar{x}$
1.17	0.71	1.41	1.24	1.7	0.93

Table 4. Table 4: Estimated parameters and population mean in the model with uncorrected selection bias. 95% Bootstrap Confidence Intervals of x ¯ ¯ 𝑥 \bar{x} is (1.38, 1.72).

$α$	$β_{1}$	$β_{2}$	$β_{3}$	$γ$	$\bar{x}$
1.80	0.71	1.41	1.25	1.54	1.46

Table 5. Table 5: Summary of the Parameters and Population Mean by Gibbs Sampler

Name	Min.	1st Qu.	Median	Mean	3rd Qu.	Max.
$α$	0.25	0.28	1.62	2.86	5.36	6.78
$β_{1}$	0.25	0.25	0.40	0.55	0.84	1.25
$β_{2}$	0.25	0.25	0.62	0.82	1.33	1.91
$β_{3}$	0.25	0.25	0.53	0.92	1.33	4.83
$γ$	0.48	0.58	0.83	1.06	1.35	2.50
$\bar{X}$	0.066	0.088	3.44	13.92	26.01	58.86

Table 6. Table 6: Summary of the Parameters and Population Mean by Random Sampler

Name	Min.	1st Qu.	Median	Mean	3rd Qu.	Max.
$α$	0.25	0.77	1.37	1.34	1.93	2.33
$β_{1}$	0.07	0.29	0.46	0.52	0.71	1.50
$β_{2}$	0.17	0.49	0.75	0.83	1.12	2.19
$β_{3}$	0.13	0.56	0.87	0.96	1.29	2.98
$γ$	0.64	1.05	1.36	1.43	1.68	3.54
$\bar{X}$	0.31	0.67	0.75	0.74	0.81	1.01

Equations105

g (x) = \frac{w ( x ) f ( x )}{\int w ( x ) f ( x ) d x} .

g (x) = \frac{w ( x ) f ( x )}{\int w ( x ) f ( x ) d x} .

g (x) = \frac{x f ( x )}{μ},

g (x) = \frac{x f ( x )}{μ},

f (x ∣ α, β, γ) = \frac{γ x ^{γ α - 1}}{β ^{γ α} Γ ( α )} exp [- (\frac{x}{β})^{γ}], x > 0,

f (x ∣ α, β, γ) = \frac{γ x ^{γ α - 1}}{β ^{γ α} Γ ( α )} exp [- (\frac{x}{β})^{γ}], x > 0,

E (x) = β \frac{Γ ( α + \frac{1}{γ} )}{Γ ( α )}, and V a r (x) = β^{2} \frac{Γ ( α + \frac{2}{γ} )}{Γ ( α )} - (\frac{Γ ( α + \frac{1}{γ} )}{Γ ( α )})^{2},

E (x) = β \frac{Γ ( α + \frac{1}{γ} )}{Γ ( α )}, and V a r (x) = β^{2} \frac{Γ ( α + \frac{2}{γ} )}{Γ ( α )} - (\frac{Γ ( α + \frac{1}{γ} )}{Γ ( α )})^{2},

g (x ∣ I = 1)

g (x ∣ I = 1)

= \frac{C x f ( x )}{\int C x f ( x ) d x} = \frac{x f ( x )}{E ( X )} .

g (x ∣ α, β, γ)

g (x ∣ α, β, γ)

= \frac{γ x ^{γ α}}{β ^{γ α + 1} Γ ( α + \frac{1}{γ} )} exp [- (\frac{x}{β})^{γ}], x > 0.

E (x) = \frac{β Γ ( α + \frac{2}{γ} )}{Γ ( α + \frac{1}{γ} )}, and V a r (x) = \frac{β ^{2} Γ ( α + \frac{3}{γ} )}{Γ ( α + \frac{1}{γ} )} - (\frac{β Γ ( α + \frac{2}{γ} )}{Γ ( α + \frac{1}{γ} )})^{2} .

E (x) = \frac{β Γ ( α + \frac{2}{γ} )}{Γ ( α + \frac{1}{γ} )}, and V a r (x) = \frac{β ^{2} Γ ( α + \frac{3}{γ} )}{Γ ( α + \frac{1}{γ} )} - (\frac{β Γ ( α + \frac{2}{γ} )}{Γ ( α + \frac{1}{γ} )})^{2} .

\overset{ˉ}{X} = \frac{1}{N} i = 1 \sum N x_{i} = f \overset{x}{ˉ}_{s} + (1 - f) \overset{ˉ}{X}_{n s},

\overset{ˉ}{X} = \frac{1}{N} i = 1 \sum N x_{i} = f \overset{x}{ˉ}_{s} + (1 - f) \overset{ˉ}{X}_{n s},

x_{ij} ∣ α, β_{i}, γ \sim in d GG (α, β_{i}, γ), j = 1, \dots, N_{i}, i = 1, \dots, ℓ,

x_{ij} ∣ α, β_{i}, γ \sim in d GG (α, β_{i}, γ), j = 1, \dots, N_{i}, i = 1, \dots, ℓ,

x_{ij} ∣ α, β_{i}, γ \sim in d GG (α + \frac{1}{γ}, β_{i}, γ), j = 1, \dots, n_{i}, i = 1, \dots, ℓ .

x_{ij} ∣ α, β_{i}, γ \sim in d GG (α + \frac{1}{γ}, β_{i}, γ), j = 1, \dots, n_{i}, i = 1, \dots, ℓ .

\hat{N} = i = 1 \sum n \frac{1}{π _{i}}

\hat{N} = i = 1 \sum n \frac{1}{π _{i}}

π_{i} = C x_{i}, i = 1, \dots, n,

π_{i} = C x_{i}, i = 1, \dots, n,

\hat{N} = 125 \times i = 1 \sum n \frac{1}{x _{i}},

\hat{N} = 125 \times i = 1 \sum n \frac{1}{x _{i}},

n_{i} ∣ N_{i}, μ_{o} \sim in d \mbox B in o mia l (N_{i}, μ_{o}), n_{i} = 0, \dots, N_{i}, i = 1, \dots, ℓ,

n_{i} ∣ N_{i}, μ_{o} \sim in d \mbox B in o mia l (N_{i}, μ_{o}), n_{i} = 0, \dots, N_{i}, i = 1, \dots, ℓ,

π (N_{i}) \propto \frac{1}{N _{i}}, N_{i} \geq n_{i} .

π (N_{i}) \propto \frac{1}{N _{i}}, N_{i} \geq n_{i} .

π (N_{i} ∣ n_{i}, μ_{o}) = \frac{( N _{i} - 1 )!}{( n _{i} - 1 )! ( N _{i} - n _{i} )!} μ_{0}^{n_{i}} (1 - μ_{0})^{N_{i} - n_{i}}, N_{i} \geq n_{i}, i = 1, \dots, ℓ,

π (N_{i} ∣ n_{i}, μ_{o}) = \frac{( N _{i} - 1 )!}{( n _{i} - 1 )! ( N _{i} - n _{i} )!} μ_{0}^{n_{i}} (1 - μ_{0})^{N_{i} - n_{i}}, N_{i} \geq n_{i}, i = 1, \dots, ℓ,

I_{j} ∣ x_{j} \sim B er (\frac{x _{j}}{W}) \leavevmode \mbox an d \leavevmode x_{j} \sim f (x_{j})

I_{j} ∣ x_{j} \sim B er (\frac{x _{j}}{W}) \leavevmode \mbox an d \leavevmode x_{j} \sim f (x_{j})

\Rightarrow π (I_{j}, x_{j}) \propto [\frac{x _{j}}{W} f (x_{j})]^{I_{j}} [(1 - \frac{x _{j}}{W}) f (x_{j})]^{1 - I_{j}}

\Rightarrow π (x_{j} ∣ I_{j} = 0) = \frac{( 1 - \frac{x _{j}}{W} ) f ( x _{j} )}{\int ( 1 - \frac{x _{j}}{W} ) f ( x _{j} ) d x _{j}} .

\displaystyle\pi(x_{n+1},...,x_{N}|\alpha,\mathchoice{\oalign{$\displaystyle\beta$\crcr\vbox to0.86108pt{\hbox{$\displaystyle\tilde{}$}\vss}}}{\oalign{$\textstyle\beta$\crcr\vbox to0.86108pt{\hbox{$\textstyle\tilde{}$}\vss}}}{\oalign{$\scriptstyle\beta$\crcr\vbox to0.86108pt{\hbox{$\scriptstyle\tilde{}$}\vss}}}{\oalign{$\scriptscriptstyle\beta$\crcr\vbox to0.86108pt{\hbox{$\scriptscriptstyle\tilde{}$}\vss}}},\gamma,x_{1},...,x_{n},N)

\displaystyle\pi(x_{n+1},...,x_{N}|\alpha,\mathchoice{\oalign{$\displaystyle\beta$\crcr\vbox to0.86108pt{\hbox{$\displaystyle\tilde{}$}\vss}}}{\oalign{$\textstyle\beta$\crcr\vbox to0.86108pt{\hbox{$\textstyle\tilde{}$}\vss}}}{\oalign{$\scriptstyle\beta$\crcr\vbox to0.86108pt{\hbox{$\scriptstyle\tilde{}$}\vss}}}{\oalign{$\scriptscriptstyle\beta$\crcr\vbox to0.86108pt{\hbox{$\scriptscriptstyle\tilde{}$}\vss}}},\gamma,x_{1},...,x_{n},N)

= j = n + 1 \prod N \frac{( 1 - \frac{x _{j}}{W} ) f ( x _{j} )}{\int ( 1 - \frac{x _{j}}{W} ) f ( x _{j} ) d x _{j}} = j = n + 1 \prod N [\frac{1 - \frac{x _{j}}{W}}{1 - \frac{μ}{W}}] f (x_{j}),

x_{ij} ∣ α, β_{i}, γ \sim in d GG (α + 1/ γ, β_{i}, γ), j = 1, \dots, n_{i}

x_{ij} ∣ α, β_{i}, γ \sim in d GG (α + 1/ γ, β_{i}, γ), j = 1, \dots, n_{i}

π (β_{i}) \propto \frac{1}{β _{i}}, i = 1, \dots, ℓ, π (α) \propto \frac{1}{( 1 + α ) ^{2}}, π (γ) \propto \frac{1}{( 1 + γ ) ^{2}} .

π (β_{i}) \propto \frac{1}{β _{i}}, i = 1, \dots, ℓ, π (α) \propto \frac{1}{( 1 + α ) ^{2}}, π (γ) \propto \frac{1}{( 1 + γ ) ^{2}} .

\pi(x_{ij},i=1,\ldots,\ell,j=n_{i}+1,\ldots\mid N_{i},\alpha,\mathchoice{\oalign{$\displaystyle\beta$\crcr\vbox to0.86108pt{\hbox{$\displaystyle\tilde{}$}\vss}}}{\oalign{$\textstyle\beta$\crcr\vbox to0.86108pt{\hbox{$\textstyle\tilde{}$}\vss}}}{\oalign{$\scriptstyle\beta$\crcr\vbox to0.86108pt{\hbox{$\scriptstyle\tilde{}$}\vss}}}{\oalign{$\scriptscriptstyle\beta$\crcr\vbox to0.86108pt{\hbox{$\scriptscriptstyle\tilde{}$}\vss}}},\gamma)=\prod_{i=1}^{\ell}\prod_{j=n_{i}+1}^{N_{i}}[\frac{1-\frac{x_{ij}}{W}}{1-\frac{\mu}{W}}]f(x_{ij}\mid\alpha,\mathchoice{\oalign{$\displaystyle\beta$\crcr\vbox to0.86108pt{\hbox{$\displaystyle\tilde{}$}\vss}}}{\oalign{$\textstyle\beta$\crcr\vbox to0.86108pt{\hbox{$\textstyle\tilde{}$}\vss}}}{\oalign{$\scriptstyle\beta$\crcr\vbox to0.86108pt{\hbox{$\scriptstyle\tilde{}$}\vss}}}{\oalign{$\scriptscriptstyle\beta$\crcr\vbox to0.86108pt{\hbox{$\scriptscriptstyle\tilde{}$}\vss}}},\gamma),

\pi(x_{ij},i=1,\ldots,\ell,j=n_{i}+1,\ldots\mid N_{i},\alpha,\mathchoice{\oalign{$\displaystyle\beta$\crcr\vbox to0.86108pt{\hbox{$\displaystyle\tilde{}$}\vss}}}{\oalign{$\textstyle\beta$\crcr\vbox to0.86108pt{\hbox{$\textstyle\tilde{}$}\vss}}}{\oalign{$\scriptstyle\beta$\crcr\vbox to0.86108pt{\hbox{$\scriptstyle\tilde{}$}\vss}}}{\oalign{$\scriptscriptstyle\beta$\crcr\vbox to0.86108pt{\hbox{$\scriptscriptstyle\tilde{}$}\vss}}},\gamma)=\prod_{i=1}^{\ell}\prod_{j=n_{i}+1}^{N_{i}}[\frac{1-\frac{x_{ij}}{W}}{1-\frac{\mu}{W}}]f(x_{ij}\mid\alpha,\mathchoice{\oalign{$\displaystyle\beta$\crcr\vbox to0.86108pt{\hbox{$\displaystyle\tilde{}$}\vss}}}{\oalign{$\textstyle\beta$\crcr\vbox to0.86108pt{\hbox{$\textstyle\tilde{}$}\vss}}}{\oalign{$\scriptstyle\beta$\crcr\vbox to0.86108pt{\hbox{$\scriptstyle\tilde{}$}\vss}}}{\oalign{$\scriptscriptstyle\beta$\crcr\vbox to0.86108pt{\hbox{$\scriptscriptstyle\tilde{}$}\vss}}},\gamma),

\pi(\alpha,\mathchoice{\oalign{$\displaystyle\beta$\crcr\vbox to0.86108pt{\hbox{$\displaystyle\tilde{}$}\vss}}}{\oalign{$\textstyle\beta$\crcr\vbox to0.86108pt{\hbox{$\textstyle\tilde{}$}\vss}}}{\oalign{$\scriptstyle\beta$\crcr\vbox to0.86108pt{\hbox{$\scriptstyle\tilde{}$}\vss}}}{\oalign{$\scriptscriptstyle\beta$\crcr\vbox to0.86108pt{\hbox{$\scriptscriptstyle\tilde{}$}\vss}}},\gamma\mid\mathchoice{\oalign{$\displaystyle x$\crcr\vbox to0.86108pt{\hbox{$\displaystyle\tilde{}$}\vss}}}{\oalign{$\textstyle x$\crcr\vbox to0.86108pt{\hbox{$\textstyle\tilde{}$}\vss}}}{\oalign{$\scriptstyle x$\crcr\vbox to0.86108pt{\hbox{$\scriptstyle\tilde{}$}\vss}}}{\oalign{$\scriptscriptstyle x$\crcr\vbox to0.86108pt{\hbox{$\scriptscriptstyle\tilde{}$}\vss}}}_{s})\propto\frac{\gamma^{n}\left(\prod\limits_{i=1}^{\ell}\prod\limits_{j=1}^{n_{i}}x_{ij}\right)^{\gamma\alpha}}{\left(\prod\limits_{i=1}^{\ell}\beta_{i}^{n_{i}}\right)^{\gamma\alpha+1}\left[\Gamma(\alpha+\frac{1}{\gamma})\right]^{n}}\exp\left[-\sum\limits_{i=1}^{\ell}\sum\limits_{j=1}^{n_{i}}{\left(\frac{x_{ij}}{\beta_{i}}\right)}^{\gamma}\right]\frac{1}{(1+\alpha)^{2}}\frac{1}{(1+\gamma)^{2}}\prod\limits_{i=1}^{\ell}\frac{1}{\beta_{i}},

\pi(\alpha,\mathchoice{\oalign{$\displaystyle\beta$\crcr\vbox to0.86108pt{\hbox{$\displaystyle\tilde{}$}\vss}}}{\oalign{$\textstyle\beta$\crcr\vbox to0.86108pt{\hbox{$\textstyle\tilde{}$}\vss}}}{\oalign{$\scriptstyle\beta$\crcr\vbox to0.86108pt{\hbox{$\scriptstyle\tilde{}$}\vss}}}{\oalign{$\scriptscriptstyle\beta$\crcr\vbox to0.86108pt{\hbox{$\scriptscriptstyle\tilde{}$}\vss}}},\gamma\mid\mathchoice{\oalign{$\displaystyle x$\crcr\vbox to0.86108pt{\hbox{$\displaystyle\tilde{}$}\vss}}}{\oalign{$\textstyle x$\crcr\vbox to0.86108pt{\hbox{$\textstyle\tilde{}$}\vss}}}{\oalign{$\scriptstyle x$\crcr\vbox to0.86108pt{\hbox{$\scriptstyle\tilde{}$}\vss}}}{\oalign{$\scriptscriptstyle x$\crcr\vbox to0.86108pt{\hbox{$\scriptscriptstyle\tilde{}$}\vss}}}_{s})\propto\frac{\gamma^{n}\left(\prod\limits_{i=1}^{\ell}\prod\limits_{j=1}^{n_{i}}x_{ij}\right)^{\gamma\alpha}}{\left(\prod\limits_{i=1}^{\ell}\beta_{i}^{n_{i}}\right)^{\gamma\alpha+1}\left[\Gamma(\alpha+\frac{1}{\gamma})\right]^{n}}\exp\left[-\sum\limits_{i=1}^{\ell}\sum\limits_{j=1}^{n_{i}}{\left(\frac{x_{ij}}{\beta_{i}}\right)}^{\gamma}\right]\frac{1}{(1+\alpha)^{2}}\frac{1}{(1+\gamma)^{2}}\prod\limits_{i=1}^{\ell}\frac{1}{\beta_{i}},

\pi(\alpha,\mathchoice{\oalign{$\displaystyle\phi$\crcr\vbox to0.86108pt{\hbox{$\displaystyle\tilde{}$}\vss}}}{\oalign{$\textstyle\phi$\crcr\vbox to0.86108pt{\hbox{$\textstyle\tilde{}$}\vss}}}{\oalign{$\scriptstyle\phi$\crcr\vbox to0.86108pt{\hbox{$\scriptstyle\tilde{}$}\vss}}}{\oalign{$\scriptscriptstyle\phi$\crcr\vbox to0.86108pt{\hbox{$\scriptscriptstyle\tilde{}$}\vss}}},\gamma\mid\mathchoice{\oalign{$\displaystyle x$\crcr\vbox to0.86108pt{\hbox{$\displaystyle\tilde{}$}\vss}}}{\oalign{$\textstyle x$\crcr\vbox to0.86108pt{\hbox{$\textstyle\tilde{}$}\vss}}}{\oalign{$\scriptstyle x$\crcr\vbox to0.86108pt{\hbox{$\scriptstyle\tilde{}$}\vss}}}{\oalign{$\scriptscriptstyle x$\crcr\vbox to0.86108pt{\hbox{$\scriptscriptstyle\tilde{}$}\vss}}}_{s})\propto\frac{1}{(1+\alpha)^{2}}\frac{1}{(1+\gamma)^{2}}

\pi(\alpha,\mathchoice{\oalign{$\displaystyle\phi$\crcr\vbox to0.86108pt{\hbox{$\displaystyle\tilde{}$}\vss}}}{\oalign{$\textstyle\phi$\crcr\vbox to0.86108pt{\hbox{$\textstyle\tilde{}$}\vss}}}{\oalign{$\scriptstyle\phi$\crcr\vbox to0.86108pt{\hbox{$\scriptstyle\tilde{}$}\vss}}}{\oalign{$\scriptscriptstyle\phi$\crcr\vbox to0.86108pt{\hbox{$\scriptscriptstyle\tilde{}$}\vss}}},\gamma\mid\mathchoice{\oalign{$\displaystyle x$\crcr\vbox to0.86108pt{\hbox{$\displaystyle\tilde{}$}\vss}}}{\oalign{$\textstyle x$\crcr\vbox to0.86108pt{\hbox{$\textstyle\tilde{}$}\vss}}}{\oalign{$\scriptstyle x$\crcr\vbox to0.86108pt{\hbox{$\scriptstyle\tilde{}$}\vss}}}{\oalign{$\scriptscriptstyle x$\crcr\vbox to0.86108pt{\hbox{$\scriptscriptstyle\tilde{}$}\vss}}}_{s})\propto\frac{1}{(1+\alpha)^{2}}\frac{1}{(1+\gamma)^{2}}

\times \frac{γ ^{n - ℓ} ( i = 1 \prod ℓ j = 1 \prod n _{i} x _{ij}^{γ α} ) \prod _{i = 1}^{ℓ} ϕ _{i}^{n_{i} (α + \frac{1}{γ}) - 1}}{[ Γ ( α + \frac{1}{γ} ) ] ^{n}} exp [- i = 1 \sum ℓ ϕ_{i} j = 1 \sum n_{i} x_{ij}^{γ}] .

\times \frac{γ ^{n - ℓ} ( i = 1 \prod ℓ j = 1 \prod n _{i} x _{ij}^{γ α} ) \prod _{i = 1}^{ℓ} ϕ _{i}^{n_{i} (α + \frac{1}{γ}) - 1}}{[ Γ ( α + \frac{1}{γ} ) ] ^{n}} exp [- i = 1 \sum ℓ ϕ_{i} j = 1 \sum n_{i} x_{ij}^{γ}] .

\phi_{i}\mid\alpha,\gamma,\mathchoice{\oalign{$\displaystyle x$\crcr\vbox to0.86108pt{\hbox{$\displaystyle\tilde{}$}\vss}}}{\oalign{$\textstyle x$\crcr\vbox to0.86108pt{\hbox{$\textstyle\tilde{}$}\vss}}}{\oalign{$\scriptstyle x$\crcr\vbox to0.86108pt{\hbox{$\scriptstyle\tilde{}$}\vss}}}{\oalign{$\scriptscriptstyle x$\crcr\vbox to0.86108pt{\hbox{$\scriptscriptstyle\tilde{}$}\vss}}}_{s}\stackrel{{\scriptstyle ind}}{{\sim}}\mbox{Gamma}\{n_{i}(\alpha+\frac{1}{\gamma}),\sum_{j=1}^{n_{i}}x_{ij}^{\gamma}\},i=1,\ldots,\ell.

\phi_{i}\mid\alpha,\gamma,\mathchoice{\oalign{$\displaystyle x$\crcr\vbox to0.86108pt{\hbox{$\displaystyle\tilde{}$}\vss}}}{\oalign{$\textstyle x$\crcr\vbox to0.86108pt{\hbox{$\textstyle\tilde{}$}\vss}}}{\oalign{$\scriptstyle x$\crcr\vbox to0.86108pt{\hbox{$\scriptstyle\tilde{}$}\vss}}}{\oalign{$\scriptscriptstyle x$\crcr\vbox to0.86108pt{\hbox{$\scriptscriptstyle\tilde{}$}\vss}}}_{s}\stackrel{{\scriptstyle ind}}{{\sim}}\mbox{Gamma}\{n_{i}(\alpha+\frac{1}{\gamma}),\sum_{j=1}^{n_{i}}x_{ij}^{\gamma}\},i=1,\ldots,\ell.

\pi(\alpha,\gamma\mid\mathchoice{\oalign{$\displaystyle x$\crcr\vbox to0.86108pt{\hbox{$\displaystyle\tilde{}$}\vss}}}{\oalign{$\textstyle x$\crcr\vbox to0.86108pt{\hbox{$\textstyle\tilde{}$}\vss}}}{\oalign{$\scriptstyle x$\crcr\vbox to0.86108pt{\hbox{$\scriptstyle\tilde{}$}\vss}}}{\oalign{$\scriptscriptstyle x$\crcr\vbox to0.86108pt{\hbox{$\scriptscriptstyle\tilde{}$}\vss}}}_{s})\propto\prod_{i=1}^{\ell}\left\{\gamma^{n_{i}-1}\frac{\left(\prod\limits_{j=1}^{n_{i}}x_{ij}^{\gamma\alpha}\right)}{\left(\sum\limits_{j=1}^{n_{i}}x_{ij}^{\gamma}\right)^{n_{i}(\alpha+\frac{1}{\gamma})}}\frac{\Gamma\{n_{i}\left(\alpha+\frac{1}{\gamma}\right)\}}{\left(\Gamma(\alpha+\frac{1}{\gamma})\right)^{n_{i}}}\right\}\frac{1}{(1+\alpha)^{2}}\frac{1}{(1+\gamma)^{2}}.

\pi(\alpha,\gamma\mid\mathchoice{\oalign{$\displaystyle x$\crcr\vbox to0.86108pt{\hbox{$\displaystyle\tilde{}$}\vss}}}{\oalign{$\textstyle x$\crcr\vbox to0.86108pt{\hbox{$\textstyle\tilde{}$}\vss}}}{\oalign{$\scriptstyle x$\crcr\vbox to0.86108pt{\hbox{$\scriptstyle\tilde{}$}\vss}}}{\oalign{$\scriptscriptstyle x$\crcr\vbox to0.86108pt{\hbox{$\scriptscriptstyle\tilde{}$}\vss}}}_{s})\propto\prod_{i=1}^{\ell}\left\{\gamma^{n_{i}-1}\frac{\left(\prod\limits_{j=1}^{n_{i}}x_{ij}^{\gamma\alpha}\right)}{\left(\sum\limits_{j=1}^{n_{i}}x_{ij}^{\gamma}\right)^{n_{i}(\alpha+\frac{1}{\gamma})}}\frac{\Gamma\{n_{i}\left(\alpha+\frac{1}{\gamma}\right)\}}{\left(\Gamma(\alpha+\frac{1}{\gamma})\right)^{n_{i}}}\right\}\frac{1}{(1+\alpha)^{2}}\frac{1}{(1+\gamma)^{2}}.

\pi(\alpha,\gamma\mid\mathchoice{\oalign{$\displaystyle x$\crcr\vbox to0.86108pt{\hbox{$\displaystyle\tilde{}$}\vss}}}{\oalign{$\textstyle x$\crcr\vbox to0.86108pt{\hbox{$\textstyle\tilde{}$}\vss}}}{\oalign{$\scriptstyle x$\crcr\vbox to0.86108pt{\hbox{$\scriptstyle\tilde{}$}\vss}}}{\oalign{$\scriptscriptstyle x$\crcr\vbox to0.86108pt{\hbox{$\scriptscriptstyle\tilde{}$}\vss}}}_{s})\propto\prod_{i=1}^{\ell}\left\{\left(\frac{g_{i}}{a_{i}}\right)^{n_{i}\alpha}\left(\frac{\gamma^{n_{i}-1}}{a_{i}^{n_{i}/\gamma}}\right)\frac{\Gamma\{n_{i}\left(\alpha+\frac{1}{\gamma}\right)\}}{n_{i}^{n_{i}(\alpha+\frac{1}{\gamma})}\left(\Gamma(\alpha+\frac{1}{\gamma})\right)^{n_{i}}}\right\}\frac{1}{(1+\alpha)^{2}}\frac{1}{(1+\gamma)^{2}}.

\pi(\alpha,\gamma\mid\mathchoice{\oalign{$\displaystyle x$\crcr\vbox to0.86108pt{\hbox{$\displaystyle\tilde{}$}\vss}}}{\oalign{$\textstyle x$\crcr\vbox to0.86108pt{\hbox{$\textstyle\tilde{}$}\vss}}}{\oalign{$\scriptstyle x$\crcr\vbox to0.86108pt{\hbox{$\scriptstyle\tilde{}$}\vss}}}{\oalign{$\scriptscriptstyle x$\crcr\vbox to0.86108pt{\hbox{$\scriptscriptstyle\tilde{}$}\vss}}}_{s})\propto\prod_{i=1}^{\ell}\left\{\left(\frac{g_{i}}{a_{i}}\right)^{n_{i}\alpha}\left(\frac{\gamma^{n_{i}-1}}{a_{i}^{n_{i}/\gamma}}\right)\frac{\Gamma\{n_{i}\left(\alpha+\frac{1}{\gamma}\right)\}}{n_{i}^{n_{i}(\alpha+\frac{1}{\gamma})}\left(\Gamma(\alpha+\frac{1}{\gamma})\right)^{n_{i}}}\right\}\frac{1}{(1+\alpha)^{2}}\frac{1}{(1+\gamma)^{2}}.

\pi(\mathchoice{\oalign{$\displaystyle\phi$\crcr\vbox to0.86108pt{\hbox{$\displaystyle\tilde{}$}\vss}}}{\oalign{$\textstyle\phi$\crcr\vbox to0.86108pt{\hbox{$\textstyle\tilde{}$}\vss}}}{\oalign{$\scriptstyle\phi$\crcr\vbox to0.86108pt{\hbox{$\scriptstyle\tilde{}$}\vss}}}{\oalign{$\scriptscriptstyle\phi$\crcr\vbox to0.86108pt{\hbox{$\scriptscriptstyle\tilde{}$}\vss}}},\alpha,\gamma\mid\mathchoice{\oalign{$\displaystyle x$\crcr\vbox to0.86108pt{\hbox{$\displaystyle\tilde{}$}\vss}}}{\oalign{$\textstyle x$\crcr\vbox to0.86108pt{\hbox{$\textstyle\tilde{}$}\vss}}}{\oalign{$\scriptstyle x$\crcr\vbox to0.86108pt{\hbox{$\scriptstyle\tilde{}$}\vss}}}{\oalign{$\scriptscriptstyle x$\crcr\vbox to0.86108pt{\hbox{$\scriptscriptstyle\tilde{}$}\vss}}}_{s})=\pi(\mathchoice{\oalign{$\displaystyle\phi$\crcr\vbox to0.86108pt{\hbox{$\displaystyle\tilde{}$}\vss}}}{\oalign{$\textstyle\phi$\crcr\vbox to0.86108pt{\hbox{$\textstyle\tilde{}$}\vss}}}{\oalign{$\scriptstyle\phi$\crcr\vbox to0.86108pt{\hbox{$\scriptstyle\tilde{}$}\vss}}}{\oalign{$\scriptscriptstyle\phi$\crcr\vbox to0.86108pt{\hbox{$\scriptscriptstyle\tilde{}$}\vss}}}\mid\alpha,\gamma,\mathchoice{\oalign{$\displaystyle x$\crcr\vbox to0.86108pt{\hbox{$\displaystyle\tilde{}$}\vss}}}{\oalign{$\textstyle x$\crcr\vbox to0.86108pt{\hbox{$\textstyle\tilde{}$}\vss}}}{\oalign{$\scriptstyle x$\crcr\vbox to0.86108pt{\hbox{$\scriptstyle\tilde{}$}\vss}}}{\oalign{$\scriptscriptstyle x$\crcr\vbox to0.86108pt{\hbox{$\scriptscriptstyle\tilde{}$}\vss}}}_{s})\pi(\alpha,\gamma\mid\mathchoice{\oalign{$\displaystyle x$\crcr\vbox to0.86108pt{\hbox{$\displaystyle\tilde{}$}\vss}}}{\oalign{$\textstyle x$\crcr\vbox to0.86108pt{\hbox{$\textstyle\tilde{}$}\vss}}}{\oalign{$\scriptstyle x$\crcr\vbox to0.86108pt{\hbox{$\scriptstyle\tilde{}$}\vss}}}{\oalign{$\scriptscriptstyle x$\crcr\vbox to0.86108pt{\hbox{$\scriptscriptstyle\tilde{}$}\vss}}}_{s}),

Peer Reviews

No public reviews on file for this paper yet. If you reviewed it on a platform where reviews are public (OpenReview, ICLR, NeurIPS, ICML), you can paste yours below so the community can read it here.

Videos

No videos yet. Explain this paper in a talk, walkthrough, or lecture? Add one.

Taxonomy

TopicsBayesian Methods and Mixture Models · Ecology and Vegetation Dynamics Studies · Forest ecology and management

Full text

**Bayesian Inference of a Finite Population Mean

Under Length-Biased Sampling**

Zhiqing Xu, Balgobin Nandram

Department of Mathematical Sciences, Worcester Polytechnic Institute

100 Institute Road, Worcester, MA 01609

( [email protected], [email protected])

Binod Manandhar

Department of Mathematics, University of Houston

([email protected])

November 17, 2018

We present a robust Bayesian method to analyze forestry data when samples are selected with probability proportional to length from a finite population of unknown size. Specifically, we use Bayesian predictive inference to estimate the finite population mean of shrub widths in a limestone quarry dominated by re-growth of mountain mahogany. The data on shrub widths are collected using transect sampling and it is assumed that the probability that a shrub is selected is proportional to its width; this is length-biased sampling. In this type of sampling, the population size is also unknown and this creates an additional challenge. The quantity of interest is average finite population shrub width and the total shrub area of the quarry can be estimated. Our method is assisted by using the three-parameter generalized gamma distribution, thereby robustifying our procedure against a possible model failure. Using conditional predictive ordinates, we show that the model, which accommodates length bias, performs better than the model that does not. In the Bayesian computation, we overcome a technical problem associated with Gibbs sampling by using a random sampler.

Keywords: Conditional predictive ordinate, Generalized gamma distribution, Gibbs sampling, Weighted distribution, Random sampling, Robustness, Transect sampling.

1. Introduction

Unequal probability sampling method was first suggested by Hansen and Hurwitz (1943) who demonstrated that the use of unequal selection probabilities frequently allowed more efficient estimators of the population total than did equal probability sampling. The sampling procedure Hansen and Hurwitz (1943) proposed was length-biased sampling. It occurs when the sample selection probabilities are correlated with the values of a study variable, e.g., size variable. This problem falls under the general umbrella of selection bias problems in survey sampling.

Line intercept sampling is a length-biased method used to study certain quantitative characteristics of objects in a region of interest. In general, objects may be of any shape and size and may possess an arbitrary spatial distribution. For example, these objects may be shrubs or patches of vegetation in a field, or the projection of logs on the forest floor. The idea of line intercept sampling is to use lines (transects) as sampling units and measuring features of the objects (e.g., widths of shrubs) that crossed by them. A length-biased sampling method producing samples from a weighted distribution. With the underlying distribution of population, one can estimate the attributes of the population by converting the weighted samples to random samples (surrogate samples).

For the estimation of a finite population quantity, the problem is more complex than for a superpopulation parameter because if there is a bias which tends to make the sampled values large, the nonsampled values would tend to be small. Such an adjustment is difficult to carry out. Generally, it has been assumed that the sample size is much smaller than the population size, and this eliminates the finite population estimation problem. Recently Nandram, Bhatta and Bhadra (2013) proposed a Bayesian non-ignorable selection model to accommodate a selection mechanism for binary data; see also Nandram (2007).

There are several approaches to address the selection bias problem. One approach incorporates the nonsampled selection probabilities in a model. This approach is computer-intensive because the nonsampled part of the population is much larger than the sample; e.g., see Nandram, Choi, Shen and Burgos (2006), Nandram and Choi (2010) and Choi, Nandram and Kim (2017). The second approach involves two models, one for the sample, called the survey model, and the other for the population, called the census model. This approach is sometimes called the surrogate sampling approach; e.g., see Nandram (2007) and Nandram, Bhatta and Bhadra (2013). The surrogate sampling approach obtains a surrogate random sample from the census model, and then prediction is done via the census model. The third approach is based on finite population sampling, in which a sample distribution and a sample complement distribution are both constructed, see Sverchkov and Pfeffermann (2004). It is convenient to use this approach for line intersect sampling. Sverchkov and Pfeffermann (2004) developed design consistent predictors for the finite population total. Essentially they define the distributions of the sampled values and the nonsampled values as two separate weighted distributions of the census distribution (see Patil and Rao 1978). Yet another approach is based on quasi-likelihood (Chambers and Skinner 2003) which is difficult to perform in a Bayesian paradigm because the normalization constant is hard to evaluate (typically a complicated function of the model parameters).

Length-biased distributions are a special case of the more general form known as weighted distributions. First introduced by Fisher (1934) to model ascertainment bias, weighted distributions were later formalized in a unifying theory by Rao (1965); see also the celebrated paper of Patil and Rao (1978). Briefly, if the random variable $x$ has a probability density function (pdf) of $f(x)$ , and a non-negative weight function of $w(x)$ , then the corresponding weighted density function is

[TABLE]

A special case is when the weight function is $w(x)=x$ . Such a distribution is known as a length-biased distribution and is given by

[TABLE]

where $\mu=\int xf(x)dx\neq 0$ , subject to existence, is the expectation of $x$ . (In Bayesian statistics we do not use upper and lower cases to differentiate random variables and fixed quantities.)

Various works have been done to characterize relationships between the original distribution and the length-biased distribution. Muttlak and McDonald (1990) suggested using a ranked set sampling procedure to estimate the population size and population mean.

In this paper, we use a three-parameter generalized gamma distribution as the original distribution to model the widths of shrubs sampled by line intercept method. Line intercept method has been found in widespread applications when estimating particle density, coverage and yields. For example, Lucas and Seber (1977) and Eberhardt (1978) derived unbiased estimators for density and percentage cover for any spatial distribution and randomly located transects. McDonald (1980) showed that the Lucas and Seber (1977) estimators for density and percentage cover are unbiased for a simple random sample of unequal length transects. Shrubs can be collected from either randomly located or systematically located transects (Butler and McDonald, 1983). It is evident that shrubs with larger widths have higher probabilities of selection.

In Section 2, we provide the background of our study. This includes data description and a introduction of the three-parameter generalized gamma distribution, which allows us to robustify our Bayesian method to accommodate the length bias. Section 3 contains the model and the results. The procedure involves the following steps. First, we derive the population size distribution as well as the sample complement distribution by Bayes’ theorem. Next, we propose a random sampling method to generate random parameters from their joint posterior distribution. Then, using each set of parameter values, we obtain a set of samples and the corresponding complement samples (Sverchkov and Pfeffermann, 2004). Finally, one sample of the finite population mean can be obtained by taking the average of the pooled samples. The goodness of fit is checked by utilizing conditional predictive ordinates, with results compared both for the models with and without the length bias. The conclusion is made in section 4.

2. Data and Robustness

2.1. Description of the Data

The data we use were collected using the line intercept sampling method (Muttlak and McDonald 1990). The study was conducted in a limestone quarry dominated by regrowth of mountain mahogany. The study area was defined by the area east of the baseline and within the walls of the quarry, where the baseline was established approximately parallel to the fissures; see Figure 1 in Appendix A. By dividing the baseline into three equal parts, three systematically placed transects were established. To ensure uniform coverage over the study area, two independent replications, each with 3 transects were selected (see Figure 1). One quantity of interest is the mean width of the shrubs in the quarry. So the variable we study is the width of the projection of the shrub encountered by transects onto the baseline (for illustration see Figure 2).

We use the data from both replications, as shown in Appendix A Tables 1 and 2. The numbers of shrubs counted in three transects are respectively $18,22,6$ and in the two transects of Replicate 2 are respectively $32,11$ (one transect of Replicate 2 has no data). Looking at the box plots of these two replications (Figure 3), we notice the clear differences in the distributions among the three transects in Replication 1; whereas in Replication 2, the distributions are close. Therefore, when making inferences using Replication 1, we regard the data from these three transects as from three different strata (they are actually so), and distinguish them in our modeling.

One complication is that we do not know the number of shrubs in the entire quarry. As we intend to use the Bayesian approach, the data from Replicate 2 is used to construct a prior distribution of the finite population size and the data from Replicate 1 is to provide inference for the population mean.

2.2. Generalized Gamma Distribution

The generalized gamma distribution (GG) was first defined by Stacy (1962) and it encompasses various subfamilies including the Weibull distribution, the generalized normal distributions and the lognormal distribution as a limit. Khodabin and Ahmadabadi (2010) provided details of the subfamilies of generalized gamma distribution. Because of the flexibility of generalized gamma distribution, we use it as the underlying population distribution in our models.

Some authors have advocated the use of simpler models because of estimation difficulties caused by the complexity of GG parameter structure. For example, Parr and Webster (1965), Hager and Bain (1971), and Lawless (1980) have considered maximum likelihood estimation in the three-parameter generalized gamma distribution. They reported problems with iterative solution of the nonlinear equations implied by the maximum likelihood method. They remarked that maximum likelihood estimators might not exist unless the sample size exceeds 400. (Our sample sizes are much smaller so we need to be careful.) In our paper, we perform Bayesian analyzes of generalized gamma distribution to overcome this issue.

The probability density of generalized gamma distribution is given by

[TABLE]

where $\alpha$ $\beta$ , $\gamma$ are all positive. It is worth noting that the mean and variance of $x$ are given by

[TABLE]

respectively. We write $x\sim GG(\alpha,\beta,\gamma)$ to denote a random variable with pdf, $f(x|\alpha,\beta,\gamma)$ defined by (2.1), which we call an unweighted generalized gamma distribution.

Note that when $\gamma=1$ , we get the standard gamma distribution, and by making $\gamma$ differ from $1$ many distributions are accommodated, thereby increasing the flexibility of the gamma distribution. It is in this sense we robustify our procedures.

The length-biased distribution of sample x is $g(x)=\displaystyle\frac{xf(x)}{E(x)},$ where $E(x)$ is the expectation of $x$ from the unweighted density function $f(x)$ . This can be easily derived as follows. Let $I$ denote the indicator variable, i.e., $I=1$ if the unit is selected and $I=0$ if the unit is not selected. Under length-biased sampling the probability that the unit has been selected given the value $x$ is $f(I=1|x)=Cx$ , where $C$ is a constant. By Bayes’ theorem, the sample pdf $g(x)$ is,

[TABLE]

For convenience we will write $g(x)$ for $g(x|I=1)$ .

Therefore, by using $GG(\alpha,\beta,\gamma)$ as the population distribution, the sample distribution is

[TABLE]

Note that $g(x)$ is also a generalized gamma distribution with parameters $\alpha_{g}=\alpha+\frac{1}{\gamma},\beta_{g}=\beta$ , and $\gamma_{g}=\gamma$ , denote by $x\sim\mbox{GG}(\alpha+\frac{1}{\gamma},\beta_{,}\gamma)$ , with mean and variance adjusted to

[TABLE]

We call it the weighted generalized gamma distribution.

3. Bayesian Methodology

In this section, we derive the population size distribution, the sample-complement distribution, as well as the posterior distribution for the parameters.

Denote $\ell$ as the number of transects, $N_{i}$ as the total number of shrubs in the $i^{th}$ transect, $N=\sum_{i=1}^{\ell}N_{i}$ . Note that all the $N_{i}$ and $N$ are unknown. Prior information about $N$ is needed to carry out a full Bayesian analysis. Denote $n_{i}$ as the number of shrubs from the $i^{th}$ transect, $n=\sum_{i=1}^{\ell}n_{i}$ is the number of samples. Let $x_{1},\ldots,x_{n}$ be the width of the sampled shrubs and $x_{n+1},\ldots,x_{N}$ be the widths for the nonsampled ones, which are to be predicted. The quantity of interest is

[TABLE]

where $f=\frac{n}{N}$ is the sample fraction and $\bar{x}_{s}=\frac{1}{n}\sum_{i=1}^{n}x_{i}$ and $\bar{X}_{ns}=\frac{1}{N-n}\sum_{i=n+1}^{N}x_{i}$ are respectively the sample and nonsample means. A posteriori inference is required for $\bar{X}_{ns}$ . It is worth noting that $\bar{X}_{ns}$ is not a sufficient statistic and cannot be derived directly from the sample. Therefore, one needs to draw $x_{n+1},\ldots,x_{N}$ to predict $\bar{X}_{ns}$ .

In many studies the population size $N$ is unknown and it must be estimated before inference can be made about $\bar{X}$ . Our application is no exception, however, this is easy to address as we have two sets of replicated samples. The second replicate (samples from the two transects are similarly distributed) can be used to construct a prior for $N$ . The first replicate (three transects need to be treated as three strata) is used to estimate the population mean shrub width. In this way “using the data twice” is avoided.

We assume that the population distributions for different strata are GG,

[TABLE]

accommodated the length bias, the sample distribution,

[TABLE]

The remaining problem is to find the distribution of the nonsampled values, $x_{ij},i=1,\ldots,\ell,j=n_{i}+1,\ldots,N_{i}$ , the so called the sample complement distribution (Sverchkov and Pfeffermann, 2004), which we will describe later.

In section 3.1 we show how to obtain the prior distribution for $N_{i}$ . In section 3.2 we describe the sample complement distribution. In section 3.3, we combine the results of 3.1 and 3.2 to derive the full Bayesian Model. In section 3.4, we study the posterior distributions in detail.

3.1. Prior Distribution of the Finite Population Size

We first find the estimate of $N$ based on the sample size. Then, estimates of $N_{i}$ can be obtained assuming proportional allocation.

The Horvitz-Thompson unbiased estimator of $N$ is

[TABLE]

where $\pi_{i}$ is the probability that the $i^{th}$ unit is selected; see Cochran (1977). Since the line intercept sampling gives the length biased data, we are actually sampling with probability proportion to width $x$ . Thus, we have

[TABLE]

where $C$ is a constant and $C=\frac{1}{W}$ , where $W=125$ (meters) is the length of the base line. Then, the estimated value of $N$ under selection bias is

[TABLE]

Using the data from the Replication 2, we have $\hat{N}=10,061$ . (Note that Replication 1 has $\ell=2$ strata and Replication 2 has $\ell=3$ strata.) Then, using proportional allocation in Replication 1, as $n_{1}=18$ , $n_{2}=22$ , $n_{3}=6$ , we have $\hat{N}_{1}=3,937$ , $\hat{N}_{2}=4,812$ , $\hat{N}_{3}=1,312$ .

Next, we assume

[TABLE]

$\mu_{o}$ does not depend on transects because of the nature of proportional allocation method.

Using independent noninformative priors for $N_{i}$ ,

[TABLE]

We derived the posterior distributions of $N_{i}$ ,

[TABLE]

which is a negative Binomial distribution with $E(N_{i}\mid n_{i},\mu_{o})=\displaystyle\frac{n_{i}}{\mu_{o}}$ . By equating $\hat{N}_{i}$ to $E(N_{i}\mid n_{i},\mu_{o})$ , we solve for the estimated value of $\mu_{o}$ , which is $\mu_{o}=\displaystyle\frac{n_{i}}{\hat{N}_{i}}=0.0046$ .

Therefore, based on Replication 2 our data-based prior distributions of the $N_{i}$ are independently negative binomial distributions with parameters $n_{i}$ and $\mu_{o}=0.0046,i=1,\ldots,\ell$ .

3.2. Sample-Complement Distribution

Next, we need to make inference about the non-sampled values. That is, we obtain the sample complement distribution (Sverchkov and Pfeffermann, 2004) and draw samples from it. We consider a single transect first and drop the transect indicator $i$ .

Let $I_{j}=1$ if $j\in s$ and $I_{j}=0$ if $j\not\in s$ , where $s$ denotes the sample set. Then

[TABLE]

Thus, the posterior sample complement distribution given all the parameters is

[TABLE]

where $f(x)$ is GG, $\mu$ is the expectation of $x$ , which is $\mu=\frac{\beta\Gamma(\alpha+\frac{1}{\gamma})}{\Gamma(\alpha)}$ . We use the sampling importance re-sampling (SIR) algorithm to perform the sampling. The SIR algorithm is ideal because $\displaystyle\prod_{j=n+1}^{N}f(x_{j})$ is a good proposal density and samples are easy to draw.

3.3. Full Bayesian Model

For the sample data our model is

[TABLE]

and the priori for $\alpha$ , $\beta_{i},i=1,\ldots,\ell$ and $\gamma$ are

[TABLE]

Note that the priors on the $\beta_{i}$ are improper and the priors on $\alpha$ and $\gamma$ are the $f(2,2)$ distributions ( $f(2,2)$ denotes the $f$ distribution with degrees of freedom being $(2,2)$ ), which are nearly noninformative (no moments exist) but proper.

The posterior sample complement distribution when incorporating all strata is

[TABLE]

where $f$ and $\mu$ are defined in the same way as (3.3).

The joint posterior density of $\alpha,\mathchoice{\oalign{$ \displaystyle\beta $\crcr\vbox to0.86108pt{\hbox{$ \displaystyle\tilde{} $}\vss}}}{\oalign{$ \textstyle\beta $\crcr\vbox to0.86108pt{\hbox{$ \textstyle\tilde{} $}\vss}}}{\oalign{$ \scriptstyle\beta $\crcr\vbox to0.86108pt{\hbox{$ \scriptstyle\tilde{} $}\vss}}}{\oalign{$ \scriptscriptstyle\beta $\crcr\vbox to0.86108pt{\hbox{$ \scriptscriptstyle\tilde{} $}\vss}}},\gamma$ given $\mathchoice{\oalign{$ \displaystyle x $\crcr\vbox to0.86108pt{\hbox{$ \displaystyle\tilde{} $}\vss}}}{\oalign{$ \textstyle x $\crcr\vbox to0.86108pt{\hbox{$ \textstyle\tilde{} $}\vss}}}{\oalign{$ \scriptstyle x $\crcr\vbox to0.86108pt{\hbox{$ \scriptstyle\tilde{} $}\vss}}}{\oalign{$ \scriptscriptstyle x $\crcr\vbox to0.86108pt{\hbox{$ \scriptscriptstyle\tilde{} $}\vss}}}_{s}=\{x_{ij},j=1,\ldots,n_{i},i=1,\ldots,\ell\}$ is

[TABLE]

The posterior density can be simplified by transforming $\beta_{i}$ to $\phi_{i}=1/\beta_{i}^{\gamma},i=1,\ldots,\ell$ . (Note that the jacobian of the transformation must be included.) Then, the joint posterior density of $\alpha,\mathchoice{\oalign{$ \displaystyle\phi $\crcr\vbox to0.86108pt{\hbox{$ \displaystyle\tilde{} $}\vss}}}{\oalign{$ \textstyle\phi $\crcr\vbox to0.86108pt{\hbox{$ \textstyle\tilde{} $}\vss}}}{\oalign{$ \scriptstyle\phi $\crcr\vbox to0.86108pt{\hbox{$ \scriptstyle\tilde{} $}\vss}}}{\oalign{$ \scriptscriptstyle\phi $\crcr\vbox to0.86108pt{\hbox{$ \scriptscriptstyle\tilde{} $}\vss}}},\gamma$ given $\mathchoice{\oalign{$ \displaystyle x $\crcr\vbox to0.86108pt{\hbox{$ \displaystyle\tilde{} $}\vss}}}{\oalign{$ \textstyle x $\crcr\vbox to0.86108pt{\hbox{$ \textstyle\tilde{} $}\vss}}}{\oalign{$ \scriptstyle x $\crcr\vbox to0.86108pt{\hbox{$ \scriptstyle\tilde{} $}\vss}}}{\oalign{$ \scriptscriptstyle x $\crcr\vbox to0.86108pt{\hbox{$ \scriptscriptstyle\tilde{} $}\vss}}}_{s}$ is

[TABLE]

This is not a standard posterior density, however, one can fit this model using Markov chain Monte Carlo methods (i.e., to obtain sample of $\alpha,\mathchoice{\oalign{$ \displaystyle\phi $\crcr\vbox to0.86108pt{\hbox{$ \displaystyle\tilde{} $}\vss}}}{\oalign{$ \textstyle\phi $\crcr\vbox to0.86108pt{\hbox{$ \textstyle\tilde{} $}\vss}}}{\oalign{$ \scriptstyle\phi $\crcr\vbox to0.86108pt{\hbox{$ \scriptstyle\tilde{} $}\vss}}}{\oalign{$ \scriptscriptstyle\phi $\crcr\vbox to0.86108pt{\hbox{$ \scriptscriptstyle\tilde{} $}\vss}}},\gamma$ ).

3.4. Further Study of the Posterior Density

One important problem we need to worry about is the posterior propriety of $\pi(\alpha,\mathchoice{\oalign{$ \displaystyle\phi $\crcr\vbox to0.86108pt{\hbox{$ \displaystyle\tilde{} $}\vss}}}{\oalign{$ \textstyle\phi $\crcr\vbox to0.86108pt{\hbox{$ \textstyle\tilde{} $}\vss}}}{\oalign{$ \scriptstyle\phi $\crcr\vbox to0.86108pt{\hbox{$ \scriptstyle\tilde{} $}\vss}}}{\oalign{$ \scriptscriptstyle\phi $\crcr\vbox to0.86108pt{\hbox{$ \scriptscriptstyle\tilde{} $}\vss}}},\gamma\mid\mathchoice{\oalign{$ \displaystyle x $\crcr\vbox to0.86108pt{\hbox{$ \displaystyle\tilde{} $}\vss}}}{\oalign{$ \textstyle x $\crcr\vbox to0.86108pt{\hbox{$ \textstyle\tilde{} $}\vss}}}{\oalign{$ \scriptstyle x $\crcr\vbox to0.86108pt{\hbox{$ \scriptstyle\tilde{} $}\vss}}}{\oalign{$ \scriptscriptstyle x $\crcr\vbox to0.86108pt{\hbox{$ \scriptscriptstyle\tilde{} $}\vss}}}_{s})$ . First, it is easy to see that

[TABLE]

Then, integrating out the $\phi_{i}$ , we get

[TABLE]

It is convenient to let $a_{i}=\sum_{j=1}^{n_{i}}x_{ij}^{\gamma}/n_{i}$ and $g_{i}=\left(\prod_{j=1}^{n_{i}}x_{ij}^{\gamma}\right)^{1/n_{i}}$ denote the arithmetic and geometric means of the $x_{ij}^{\gamma},j=1,\ldots,n_{i},i=1,\ldots,\ell$ . Then, we have

[TABLE]

Thus, we essentially have a two-parameter posterior density. Although an overkill, we attempted to fit this model using a Gibbs sampler. There are difficulties in performing the Gibbs sampler (perhaps associated with the difficulties encountered in finding MLEs in generalized gamma distribution) because high correlations are present among the parameters and thinning is not helpful. The problem is essentially high correlations between $\alpha$ and $\gamma$ . Thus, we consider an alternative algorithm which simply uses the multiplication rule of probability. We prove the theorem below which adds credence to our Bayesian methodology.

However, since $\gamma=1$ makes the generalized gamma density a standard gamma density, it is sensible to bound $\gamma$ in an interval centered at $1$ . That is, we take $a_{o}^{-1}\leq\gamma\leq a_{o}$ ; a sensible choice is $a_{0}=10$ or so. Thus, we replace the prior on $\gamma$ by $\gamma\sim\mbox{Uniform}(a_{o}^{-1},\leavevmode\nobreak\ a_{o})$ ; the original prior is inconvenient and not helpful.

Theorem

Assuming that $a_{o}^{-1}\leq\gamma\leq a_{o}$ , the joint posterior density of $\pi(\alpha,\gamma\mid\mathchoice{\oalign{$ \displaystyle x $\crcr\vbox to0.86108pt{\hbox{$ \displaystyle\tilde{} $}\vss}}}{\oalign{$ \textstyle x $\crcr\vbox to0.86108pt{\hbox{$ \textstyle\tilde{} $}\vss}}}{\oalign{$ \scriptstyle x $\crcr\vbox to0.86108pt{\hbox{$ \scriptstyle\tilde{} $}\vss}}}{\oalign{$ \scriptscriptstyle x $\crcr\vbox to0.86108pt{\hbox{$ \scriptscriptstyle\tilde{} $}\vss}}}_{s})$ is proper.

Remark: Using the multiplication rule of probability,

[TABLE]

the theorem implies that $\pi(\mathchoice{\oalign{$ \displaystyle\phi $\crcr\vbox to0.86108pt{\hbox{$ \displaystyle\tilde{} $}\vss}}}{\oalign{$ \textstyle\phi $\crcr\vbox to0.86108pt{\hbox{$ \textstyle\tilde{} $}\vss}}}{\oalign{$ \scriptstyle\phi $\crcr\vbox to0.86108pt{\hbox{$ \scriptstyle\tilde{} $}\vss}}}{\oalign{$ \scriptscriptstyle\phi $\crcr\vbox to0.86108pt{\hbox{$ \scriptscriptstyle\tilde{} $}\vss}}},\alpha,\gamma\mid\mathchoice{\oalign{$ \displaystyle x $\crcr\vbox to0.86108pt{\hbox{$ \displaystyle\tilde{} $}\vss}}}{\oalign{$ \textstyle x $\crcr\vbox to0.86108pt{\hbox{$ \textstyle\tilde{} $}\vss}}}{\oalign{$ \scriptstyle x $\crcr\vbox to0.86108pt{\hbox{$ \scriptstyle\tilde{} $}\vss}}}{\oalign{$ \scriptscriptstyle x $\crcr\vbox to0.86108pt{\hbox{$ \scriptscriptstyle\tilde{} $}\vss}}}_{s})$ is also proper. Clearly, $\pi(\mathchoice{\oalign{$ \displaystyle\beta $\crcr\vbox to0.86108pt{\hbox{$ \displaystyle\tilde{} $}\vss}}}{\oalign{$ \textstyle\beta $\crcr\vbox to0.86108pt{\hbox{$ \textstyle\tilde{} $}\vss}}}{\oalign{$ \scriptstyle\beta $\crcr\vbox to0.86108pt{\hbox{$ \scriptstyle\tilde{} $}\vss}}}{\oalign{$ \scriptscriptstyle\beta $\crcr\vbox to0.86108pt{\hbox{$ \scriptscriptstyle\tilde{} $}\vss}}},\alpha,\gamma\mid\mathchoice{\oalign{$ \displaystyle x $\crcr\vbox to0.86108pt{\hbox{$ \displaystyle\tilde{} $}\vss}}}{\oalign{$ \textstyle x $\crcr\vbox to0.86108pt{\hbox{$ \textstyle\tilde{} $}\vss}}}{\oalign{$ \scriptstyle x $\crcr\vbox to0.86108pt{\hbox{$ \scriptstyle\tilde{} $}\vss}}}{\oalign{$ \scriptscriptstyle x $\crcr\vbox to0.86108pt{\hbox{$ \scriptscriptstyle\tilde{} $}\vss}}}_{s})$ is also proper.

Proof

We make two observations. First, using the arithmetic-geometric inequality, we have $\left(\frac{g_{i}}{a_{i}}\right)^{n_{i}\alpha}\leq 1$ . Second, using $a_{o}^{-1}\leq\gamma\leq a_{o}$ , we have $\frac{\gamma^{n_{i}-1}}{a_{i}^{n_{i}/\gamma}}$ , a function of $\gamma$ , is bounded uniformly in $\gamma$ . Therefore, we only need to show that

[TABLE]

Next, we transform $\alpha$ to $\theta=\alpha+\frac{1}{\gamma}$ , keeping $\gamma$ untransformed. Then, the integral becomes

[TABLE]

where

[TABLE]

We only need to show that

[TABLE]

is bounded uniformly in $\theta$ for any $i=1,\ldots,\ell$ . For convenience, we will drop the subscript, $i$ , momentarily, so we simply need to show that $\Delta(\theta)=\tilde{\Gamma}(n\theta)-n\tilde{\Gamma}(\theta)-n\theta\ln(n)$ , where $\tilde{\Gamma}(\cdot)$ is the logarithm of the gamma function, is uniformly bounded in $\theta$ ; see the Appendix A.

Finally, because $g^{\ast}(\theta)\leq A<\infty$ , we are left with

[TABLE]

Therefore, our claim on propriety holds.

4. Bayesian Computations and Data Analyzes

In this section, we perform Bayesian analysis of the posterior distributions of population parameters by a numerical method, called random sampler, which performs better than the Gibbs sampler . We then obtain the nonsampled values using the sampling importance resampling (SIR) algorithm. Recall that the data we use here is Replication I, which has three transects, i.e., $\ell=3$ .

4.1. Random Sampler

Since the Gibbs sampler is a Markovian updating scheme, we have shown, in our work not presented in this paper, that most of the estimated values of population mean are larger than what we expected (see results in Appendix Table 3). One of the reasons is that high correlations among these parameters makes the Gibbs sampler inefficient in the sense it may take a very large number of iterations to converge in distribution. In this section, we propose a non-Markovian algorithm, called random sampler, in order to avoid the particular issue mentioned above.

Therefore, $\alpha$ and $\gamma$ cannot be sampled directly from their unbounded parameters space. We use the transformation $\alpha^{\prime}=\frac{\alpha}{1+\alpha}$ and $\gamma^{\prime}=\frac{\gamma}{1+\gamma}$ . Then,

[TABLE]

Two-dimensional grid method can be applied to draw $\alpha^{\prime}$ and $\gamma^{\prime}$ from their joint distribution. But grid method is computationally intensive in more than one dimension. We used the multiplication rule to draw samples of $\alpha^{\prime}$ and $\gamma^{\prime}$ .

[TABLE]

To apply this rule, we fist generated a sample of ${\gamma^{\prime}}^{(1)}$ from $\Pi(\gamma^{\prime}|x_{11}\cdots,x_{3n_{3}})$ , then generated a sample of ${\alpha^{\prime}}^{(1)}$ from $\Pi(\alpha^{\prime}|{\gamma^{\prime}}^{(1)},x_{11}\cdots,x_{3n_{3}})$ . Repeating this procedure $M$ times to obtain $M$ sets of $\alpha^{\prime}(\alpha)$ and $\gamma^{\prime}(\gamma)$ . The corresponding $\phi(\beta)$ can also be obtained by sampling from $\Pi(\phi_{i}|\alpha,\phi_{k},\gamma,x_{11},\cdots,x_{3n_{3}})$ .

The term $\Pi(\alpha^{\prime}|\gamma^{\prime},x_{11}\cdots,x_{3n_{3}})$ in (4.1) is easy to derive.

[TABLE]

The term $\Pi(\gamma^{\prime}|x_{11}\cdots,x_{3n_{3}})$ in (4.1) can be derived by integrating $\Pi(\alpha^{\prime},\gamma^{\prime}|x_{11}\cdots,x_{3n_{3}})$ with respect to $\alpha^{\prime}$ . Unfortunately, it is not possible to integrate $\Pi(\alpha^{\prime},\gamma^{\prime}|x_{11}\cdots,x_{3n_{3}})$ by analytical techniques. For this reason, numerical methods have to be used. We use the 20-point Gaussian quadrature to approximate $\Pi(\gamma^{\prime}|x_{11}\cdots,x_{3n_{3}})$ .

[TABLE]

where $x_{i},i=1,\cdots,20$ are the roots of orthogonal polynomials $P_{20}(x)$ for $[-1,1]$ and $\omega_{i},i=1,\cdots,20$ are the corresponding Gauss-Legender weights, which can be created by R package ’gaussquad’.

The summary of samples drawn by Random Sampler for each parameter and the population mean are shown in Table 6 . It is shown that the population mean has the IQR of $(.67,.81)$ , with the median of $.75$ . In Fig. 4, it is shown that the posterior distribution of $\alpha$ is bimodal (pointing to the difficulty in estimating $\alpha$ ), while those of $\beta_{1}$ , $\beta_{2}$ , $\beta_{3}$ and $\gamma$ are skewed to the right. In the next section, we will perform the model checking by conditional predictive ordinate (CPO).

4.2. Model Checking by Conditional Predictive Ordinate

Comparing the predictive distribution to the observed data is generally termed a “posterior predictive check”. This type of check includes the uncertainty associated with the estimated parameters of the model. Posterior predictive checks (via the predictive distribution) involve a double-use of the data, which causes predictive performance to be overestimated. To overcome this drawback, Geisser and Eddy (1979) has proposed the leave-one-out cross-validation predictive density. This is also known as the conditional predictive ordinate or CPO (Gelfand, 1996). The CPO is a handy posterior predictive check because it may be used to identify outliers, influential observations, and for hypothesis testing across different non-nested models. The CPO expresses the posterior probability of observing the value of $x_{i}$ when the model is fitted to all data except $x_{i}$ , with a larger value implying a better fit of the model to $x_{i}$ , and very low CPO values suggest that $x_{i}$ is an outlier and an influential observation. A Monte Carlo estimate of the CPO is obtained without actually omitting $x_{i}$ from the estimation, and is provided by the harmonic mean of the likelihood for $x_{i}$ . Specifically, the CPOi is the inverse of the posterior mean of the inverse likelihood of $x_{i}$ . The Monte Carlo estimate of CPO is

[TABLE]

where $\tilde{\theta}^{(h)}\overset{iid}{\sim}\Pi(\tilde{\theta}|\tilde{x})$ for $h=1,\dots,M$ ; see Molina, Nandram and Rao (2014).

The sum of the log(CPO)’s can be an estimator for the natural logarithm of the marginal likelihood, sometimes called the log pseudo marginal likelihood (LPML)

[TABLE]

Models with larger LPMLs are better. To compare the predictive distributions (both model with length bias and model without length bias) using our length-biased sample, we calculated the LPML for both models. The likelihood of $x_{i}$ under length biased model is given by

[TABLE]

where $\beta$ is the corresponding parameter for the stratum that $x_{i}$ is from. The likelihood of $x_{i}$ under no length biased model is given by

[TABLE]

where $\beta$ is the corresponding parameter for the stratum that $x_{i}$ is from. It is shown that the $LPML$ of for model with length bias is larger (LPML = $-36.10$ ) than the one for the model without length bias (LPML = $-47.54$ ), which means the model with length bias fits our length-biased sample better.

4.3. Nonsampled Widths

We define the importance function as

[TABLE]

Then, the importance ratios are

[TABLE]

A random sample can now be obtained by re-sampling with probability proportional to the ratios.

The algorithm to obtain the nonsampled values is as follows.

•

Step 1. Obtain M sets of $(\alpha,\beta_{1},\beta_{2},\beta_{3},\gamma,)$ using the sampling methods described in the next Section.

•

Step 2. Obtain a sample of $N$ from formula (3.2).

•

Step 3. For each set of parameters, generate the vector $\tilde{x_{j}}$ where $x_{ij},i=n_{j}+1,\cdots,(N_{j}-n_{j}),\;j=1,2,3$ from the corresponding generalized gamma distribution.

•

Step 4. Compute the population mean and the importance ratio $w$ .

•

Step 5. Repeat the step 2 to step 4 $M-1$ times.

•

Step 6. Draw $\alpha M$ values of the population mean with probability proportional to $\tilde{w}$ . We choose $\alpha=.1$

5. Summary

In this paper we have presented a model for estimating population mean under length-biased sampling. We have used a weighted distribution of the three-parameter generalized gamma distribution to model the shrub widths to robustify our procedure that accommodates the length-biased sampling. Interest is on the finite population mean of shrub width in the entire quarry. In order to avoid certain technical issues associated with classical inference when using the generalized gamma distribution, we proposed a non-Markovian Bayesian numerical method, called random sampler, which performs better than Gibbs sampler when the population parameters are highly correlated. Posterior population distributions are easily estimated using this method. Conditional predictive ordinate shows that the model with length bias performs better than the model without length bias.

While accommodating transect sampling is a challenge, another important challenge in our procedure is to estimate the unknown population size. To ensure a full Bayesian procedure, we have used the data in Replication 2 (two strata) to estimate the finite population size that is usually unknown in this type of problem. The data from Replication 1 were used to estimate the average shrub width of the finite population. This estimate can, in turn, be used to give an estimate of the total shrub area assuming a standard geometry (e.g., a circle with the width being the diameter or a square with the width being length of a side).

An interesting topic for future research would be including covariates to study potential predictors. In Muttlak and McDonald (1990), in addition to the measurement of shrub widths, two more attributes of mountain mahogany, maximum height, and number of stems, were measured. Both attributes are important predictors of the average shrub width of an area’s vegetation. Semi-parametric linear regression (Chen, 2010) or generalized linear regression can be considered to measure this association. We can incorporate the covariates through a gamma type regression model. Let the covariates be $\underset{\sim}{z_{ij}},i=1,2,3,j=n_{1},n_{2},n_{3}$ and $\underset{\sim}{\phi_{i}},i=1,2,3$ . Because the mean of each stratum is linearly related to $\beta_{1},\beta_{2},\beta_{3}$ respectively, we take $\beta_{i}=\displaystyle e^{\underset{\sim}{z_{i}}^{\prime}\underset{\sim}{\phi_{i}}},i=1,2,3$ . For the shrub data, our model is

[TABLE]

A similar form can be easily written down for the length-biased sampling. Our future plan is to fit a model to accommodate the covariates.

Acknowledgement.

*The author thanks the Associate Editor, XXX, and the referee for their comments which improved the quality of the paper very much. *

**Appendix A. Uniform Boundedness of $\Delta(\theta)$

**

We need to show that

[TABLE]

is uniformly bounded in $\theta$ . We will show that $\Delta(\theta)$ asymptotes out horizontally.

First, differentiating $\Delta(\theta)$ , we have

[TABLE]

where $\psi(\cdot)$ is the diagamma function. Now, using the duplication property (Abramowitz and Stegun 1965, Ch. 6) of the digamma function, one can show that $\Delta^{\prime}(\theta)\geq 0$ . That is, $\Delta(\theta)$ is monotonically increasing in $\theta$ ; see also Figure 5.

Second, differentiating $\Delta^{\prime}(\theta)$ , we have

[TABLE]

Using a theorem (Ronning 1986) which states that $x\psi^{\prime}(x)$ decreases monotonically in $x$ , we have $\Delta^{\prime\prime}(\theta)\leq 0$ . That is, $\Delta(\theta)$ is concave and the rate of increase of $\Delta(\theta)$ decreases; see Figure 6.

Therefore, $\Delta(\theta)$ asymptotes out horizontally and $\Delta(\theta)$ must be bounded; so is its exponent.

Bibliography30

The reference list from the paper itself. Each links out to its DOI / PubMed record.

1[1]
22 Abramowitz, M. and Stegun, I. A. (1965), Handbook of Mathematics Functions , Dover Publications.
33 Butler, S. A. and Mc Donald, L. L. (1983), “ Unbiased systematic sampling plans for the line intercept method,” Journal of Range Management , 36, 463-468.
44 Chambers, R. L. and Skinner, C. J. (2003), Analysis of Survey Data , Wiley.
55 Chen, Y. Q. (2010), “Semiparametric regression in size-biased sampling, ” Biometrics , 66, 149-158.
66 Choi, S., Nandram, B. and Kim, D. (2017), “A hierarchical Bayesian model for binary data incorporating selection bias, ” Communications in Statistics: Simulation and Computation , 46 (6), 4767-4782.
77 Cochran, W. G. (1977), Sampling Techniques , Wiley, New York.
88 Eberhnrdt, L. L. (1978), “Transect methods for population studies, ” The Journal of Wildlife Management , 42, 1-31.

TL;DR

Contribution

Findings

Abstract

Peer Reviews

Videos

Taxonomy

1. Introduction

2. Data and Robustness

2.1. Description of the Data

2.2. Generalized Gamma Distribution

3. Bayesian Methodology

3.1. Prior Distribution of the Finite Population Size

3.2. Sample-Complement Distribution

3.3. Full Bayesian Model

3.4. Further Study of the Posterior Density

4. Bayesian Computations and Data Analyzes

4.1. Random Sampler

4.2. Model Checking by Conditional Predictive Ordinate

4.3. Nonsampled Widths

5. Summary

Acknowledgement**.**

Acknowledgement.