---
title: "Getting Started with SHELF"
output: rmarkdown::html_vignette
vignette: >
  %\VignetteIndexEntry{Getting Started with SHELF}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---

```{r, include = FALSE}
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)
```

# Introduction

We illustrate the process of fitting a probabilty distribution to an expert's judgements, and obtaining feedback. 

Various elicitation methods are suggested in SHELF, but they all involve either asking for probability or quantile judgements: the expert is either asked questions of the form, 'What is your probability $P(X\le10)$?' or, 'Provide a value $x$ such that you think $P(X\le x)=0.5.$'' Either way, we think of the expert as providing points on his/her cumulative distribution function (CDF). 

As an example, suppose the uncertain quantity is the percentage $X$ of voters who will vote for a particular candidate in an upcoming leadership election. 

# Eliciting individual distributions from multiple experts

Suppose we have two experts, who first consider their judgements independently. Expert A states $$P(X\le25)=0.25, P(X\le30)=0.5, P(X\le35)=0.75,$$ and
Expert B states $$P(X\le30)=0.25, P(X\le35)=0.5, P(X<50)=0.75.$$



The values of the quantiles are arranged in a matrix, with one column per expert
```{r define v}
v <- matrix(c(25, 30, 35, 30, 35, 50), nrow = 3, ncol = 2)
```
```{r show v, echo = F}
v
```

The probabilities should be arranged in matrix if they are different for the two experts (for example, if one expert had specified tertiles instead of lower and upper quartiles), otherwise they are stored in a vector

```{r define p}
p <- c(0.25, 0.5, 0.75)
```

We now fit probability distributions to these judgements with the `fitdist` command. Without specifying any limits (`lower` and `upper`), only normal and Student-$t$ distributions (with 3 degrees of freedom used as a default) will be fitted. Including a `lower` limit will result in log-normal, log Student-$t$ and gamma distributions also being fitted, and including both a `lower` and `upper` limit will result in a beta distribution fit, scaled to the interval [`lower`, `upper`].
```{r use myfit}
library(SHELF)
myfit <- fitdist(vals = v, probs = p, lower = 0, upper = 100)
```

The object `myfit` includes parameters of the fitted distributions for each expert, sums of squared errors (elicited - fitted probabilities) for each distribution and expert, and the original elicited judgements:
```{r show myfit}
names(myfit)
```
For example, the parameters of the fitted beta distributions are 
```{r show beta}
myfit$Beta
```

Inspecting `myfit$ssq` (sum of squared differences between elicited and fitted probabilities), we see that the normal distribution fits best for Expert 1, and the log Student-$t_3$ distribution fits best for Expert 2 (although the fit would probably be judged adequate for any the distributions, in this case). These best-fitting distributions are used as defaults in the `plotfit` function, and the `feedback` function when used with multiple experts

```{r show ssq}
myfit$ssq
```

We plot the fitted distributions, including a linear pool.
```{r plot-fitted-distributions, fig.height = 4, fig.width = 5, fig.pos="h", fig.align="center", fig.cap = "The two fitted distributions and an equal-weighted linear pool."}
plotfit(myfit, lp = TRUE)
```

# Eliciting a single distribution from an individual or a group of experts

Now we elicit a single a single 'consensus' distribution from the two experts.
Suppose they agree $$P(X\le 25)=0.25, P(X\le30)=0.5, P(X\le40)=0.75.$$
We fit a single set of distributions to these judgements.
```{r define single set}
v <-c(25, 30, 40)
p <-c(0.25, 0.5, 0.75)
consensus <- fitdist(vals = v, probs = p, lower = 0, upper = 100)
```

The best fitting distribution is the log Student-$t_3$ distribution, but the (scaled) beta distribution would likely be adequate at this stage, and we may choose to use the beta instead, given it has the appropriate bounds.

We plot the fitted density, and mark the 
fitted 5th and 95th percentiles
```{r plot-RIO, fig.height = 4, fig.width = 5, fig.pos="h", fig.align="center", fig.cap = "The fitted consensus distribution, with the lower and upper 5\\% tail areas shown as feedback."}
plotfit(consensus, ql = 0.05, qu = 0.95, d = "beta")
```

We can obtain values of fitted percentiles and probabilities with the `feedback` function. This will show fitted values for all the fitted distributions 
```{r show feedback}
feedback(consensus, quantiles = c(0.05, 0.95))
```
For the fitted beta distribution, we have $P(X<15.2)=0.05$, and we can also compare, for example, the fitted probability $P(X<25)=0.283$ with elicited probability $P(X<25)=0.25$.