explodecomputer

library(ggplot2)
library(reshape2)
library(dplyr)
library(latex2exp)


# Read in data
dat <- read.table("Downloads/data_for_figures.txt", he=T, sep=",")

explodecomputer

# Testing metaCCA estimation of phenotypic correlation

## Introduction

Phenotypic correlation between two traits is a function of the shared genetic effects and shared environmental effects. For example suppose two traits are influenced by a single SNP, and the effect size of the SNP is the same on both traits. The genetic correlation is 1. If they also share the same environmental factors, but the effect of the factor on the first trait is positive, and the effect on the second trait is negative, the environmental correlation is -1. The phenotypic correlation of these two traits would be 0 (assuming equal variance of genetic and environmental factors).

These simulations just show that the metaCCA method of estimating phenotypic correlations are actually only estimating the genetic correlations. It also compares the estimate using a single sample or two samples.

From the perspective of metaCCA, whose objective is to estimate the joint effects of a SNP on correlated outcomes, it doesn't matter that it is e

explodecomputer

library(gtools)
library(dplyr)

func <- function(g, s, b)
{
	l <- permutations(3, 3, v=c(g,s,b))
	return(l)
}

l <- func(5, 2, 1)

explodecomputer

# Setup bash profile

mkdir repo
cd repo
git clone [email protected]:explodecomputer/bashInit.git
cd bashInit
./init.sh


# Mount hard drive

explodecomputer

n <- 100000
nsnp1 <- 10
nsnp2 <- 10
G1 <- matrix(rbinom(n*nsnp1, 2, 0.5), n, nsnp1)
G2 <- matrix(rbinom(n*nsnp2, 2, 0.5), n, nsnp2)

G2[,1] <- G1[,1]

eff1 <- rnorm(nsnp1)
eff2 <- rnorm(nsnp2)

explodecomputer

library(readxl)
a <- read_excel("~/Downloads/CVD_traits_summary.xlsx", skip=2)
names(a)[1] <- "Phenotype"

a1 <- a[,c(1,2:4)]
a2 <- a[,c(1,5:7)]
a3 <- a[,c(1,8:10)]
a4 <- a[,c(1,11:13)]
a5 <- a[,c(1,14:16)]
a6 <- a[,c(1,17:19)]

explodecomputer

library(ggplot2)
library(ggthemes)
library(readxl)

b <- read_excel("~/Downloads/CAD-metabolites-3-way-comparison.xlsx", skip=1)
b <- b[!is.na(b$exposure), ]
names(b) <- gsub("-", "\\.", names(b))
b$obs.logOR <- log(b$obs.OR)
b$obs.SE <- NA

explodecomputer

n <- 1000
g <- rbinom(n, 2, 0.5)
y <- g + rnorm(n)
yb <- rep(0, n)
yb[y >= median(y)] <- 1

dat <- data.frame(y=y, yb=yb, g=g)
dat <- dat[order(dat$yb), ]
dat$group <- rep(1:sum(yb==0), 2)

explodecomputer

spring

explodecomputer

library(WeightedCluster)

data(mvad)
mvad.alphabet <- c("employment", "FE", "HE", "joblessness", "school", "training")
mvad.labels <- c("Employment", "Further Education", "Higher Education", "Joblessness", "School", "Training")
mvad.scodes <- c("EM", "FE", "HE", "JL", "SC", "TR")
mvadseq <- seqdef(mvad[, 17:86], alphabet = mvad.alphabet,
	states = mvad.scodes, labels = mvad.labels,
	weights = mvad$weight, xtstep = 6)