Monday, October 23, 2023

SNP effects from Single Step GBLUP with APY

 APY is a technique that allows representing the inverse of a genomic relationship matrix in a sparse format by choosing a "core" of animals (here and here). The authoritative guide to APY is Bermann et al. 2022: this paper

One of the key aspects of genomic models is the ability of estimating SNP effects and then apply them to newly genotyped animals, what is know as Indirect Predictions. Matias and I found out that it's easier than we thought as shown here.

Among many other things Bermann et al. show that one can write indirect predictions of "non-core" animals as (I use the original equation numbering)

eq. (10) : $latex \mathbf{{u}_n}= \mathbf{Z}_n \mathbf{Z}'_c (\mathbf{Z}_c \mathbf{Z}'_c)^{-1} \mathbf{Z}_c  \mathbf{{a}} + \boldsymbol{\xi} $latex

where $latex \mathbf{{a}} $latex are SNP effects and $latex \boldsymbol{\xi} $latex is an error term that does not depend on $latex \mathbf{{u}_c} $latex. 

we obtain SNP effects from eq. 21 in Bermann et al:

$latex \mathbf{\hat{a}}=k \mathbf{Z}'_c \mathbf{G}_{cc}^{-1} \mathbf{\hat{u}}_{c}  $latex

we plug that into (10) and expand:

$latex \mathbf{\hat{u}_n}= \mathbf{Z}_n \mathbf{Z}'_c (\mathbf{Z}_c \mathbf{Z}'_c)^{-1} \mathbf{Z}_c \mathbf{\hat{a}} =k \mathbf{Z}_n \mathbf{Z}'_c (k \mathbf{Z}_c \mathbf{Z}'_c)^{-1} \mathbf{Z}_c \mathbf{\hat{a}} $latex

we substitute for $latex \mathbf{\hat{a}} $latex:

$latex \mathbf{\hat{u}_n} =k \mathbf{Z}_n \mathbf{Z}'_c (k \mathbf{Z}_c \mathbf{Z}'_c)^{-1} \mathbf{Z}_c k \mathbf{Z}'_c \mathbf{G}_{cc}^{-1} \mathbf{\hat{u}}_{c}=k \mathbf{Z}_n \mathbf{Z}'_c \mathbf{G}_{cc}^{-1} \mathbf{\hat{u}}_{c} =  \mathbf{Z}_n \mathbf{\hat{a}} $latex

which is the original eq. 21. This is because $latex \hat{a} $latex  is part of the  column space of  $latex \mathbf{Z}'_c $latex, then $latex \mathcal{P}\mathbf{\hat{a}} = \mathbf{\hat{a}} $latex.

Finally, obtaining Indirect Predictions from APY is deadly simple and intuitive:

  • $latex \mathbf{\hat{a}}=k \mathbf{Z}'_c \mathbf{G}_{cc}^{-1} \mathbf{\hat{u}}_{c}  $latex
  • $latex \mathbf{\hat{u}_n} =  \mathbf{Z}_n \mathbf{\hat{a}} $latex



Thursday, May 25, 2023

compiling macs in MacBook with Ventura

 I wanted to compile the coalescent simulator macs in the Mac running Ventura with M2 chips. This was of course the beginning of a great adventure. I end up doing the following:


- install library boost using homebrew

- dig out the path for the boost library e.g. as here, which turned out to be /opt/homebrew/Cellar/boost/1.81.0_1/include

- finally, modify the makefile as follows

# compile options
CFLAGS = -Wall -g
#CFLAGS = -Wall -O3
# Add location of any library locations below with -L
LINKFLAGS =
#LINKFLAGS = -static

# compiler
CC = g++

# libraries. For a local Boost installation
# Example:
#LIB = -I /Users/garychen/software/boost_1_36_0
# Default:
#LIB = -I .
LIB = -I /opt/homebrew/Cellar/boost/1.81.0_1/include

# simulator name
SIM=macs

OBJS = simulator.o algorithm.o datastructures.o

$(SIM) : $(OBJS)
$(CC) -o $(SIM) $(OBJS) $(LINKFLAGS)

simulator.o: simulator.cpp simulator.h
$(CC) $(CFLAGS) $(LIB) -c $<

algorithm.o: algorithm.cpp simulator.h
$(CC) $(CFLAGS) $(LIB) -c $<

datastructures.o: datastructures.cpp simulator.h
$(CC) $(CFLAGS) $(LIB) -c $<

Thursday, February 2, 2023

draw many elements julia

 # a vector of frequencies

p=rand(Beta(2,2),10)

# a vector of "distributions"

a=Binomial.(2,p)

# a vector of vectors

 b=rand.(a,1)

# collapsed (row vector)

reduce(hcat,b)


# draw 5 animals at once


b=rand.(a,5)

# collapse

bb=reduce(hcat,b)


#compute frequencies

pest = mean(bb,dims=1)/2