a02_proportionality.Rmd

---
title: "Coloc: proportionality testing"
author: "Chris Wallace"
date: "`r Sys.Date()`"
output: rmarkdown::html_vignette
vignette: >
  %\VignetteIndexEntry{Coloc: proportionality testing}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---
# Proportional testing

```{r  sim, child="simdata.Rmd" }
```

## Analysis
The code below first prepares a principal component object by combining
the genotypes in the two dataset, then models the most informative
components (the minimum set required to capture 80% of the genetic
variation) in each dataset, before finally testing whether there is
colocalisation between these models.

```{r  fig=TRUE }
## run a coloc with pcs
library(coloc)
pcs <- pcs.prepare(data@df1[,-1], data@df2[,-1])
pcs.1 <- pcs.model(pcs, group=1, Y=data@df1[,1], threshold=0.8)
pcs.2 <- pcs.model(pcs, group=2, Y=data@df2[,1], threshold=0.8)
ct.pcs <- coloc.test(pcs.1,pcs.2)
plot(ct.pcs)
```

The plot shows the estimated coefficients for each principal component
modeled for traits 1 and 2 on the x and y axes, with crosses showing
the 95% confidence intervals.  The points lie close to the line through
the origin, which supports a hypothesis of colocalisation.

A little more information is stored in the `ct.pcs` object:

```{r   }
ct.pcs
str(summary(ct.pcs)) 
```

The best estimate for the coefficient of proportionality,
\(\hat{\eta}\), is 1.13, and the null hypothesis of colocalisation is
not rejected with a chisquare statistic of 5.27 based on 7 degrees of
freedom (\(n-1\) where the \(n\) is the number of components tested, and
one degree of freedom was used in estimating \(\eta\)), giving a p value
of 0.63.  The `summary()` method returns a named vector of length 4
containing this information.

If more information is needed about \(\eta\), then this is available if
the `bayes` argument is supplied:

```{r   }
ct.pcs.bayes <- coloc.test(pcs.1,pcs.2, bayes=TRUE)
plot(ct.pcs.bayes)
ci(ct.pcs.bayes) 
```


<a id="orgb852386"></a>

## Using Bayes Factors to compare specific values of \(\eta\)

It may be that specific values of \(\eta\) are of interest.  For
example, when comparing eQTLs in two tissues, or when comparing risk
of two related diseases, the value \(\eta=1\) is of particular
interest.  In proportional testing, we can use Bayes Factors to
compare the support for different values of \(\eta\).  Eg

```{r   }
## compare individual values of eta
ct.pcs <- coloc.test(pcs.1,pcs.2, bayes.factor=c(-1,0,1))
bf(ct.pcs)

```