# 1) otcovský model/ sire model

y <- matrix(c(9,12,11,6,7,14),6,1)
y

X <- matrix(c(1,0,1,1,1,0,
              0,1,0,0,0,1),6,2)
X

Z <- matrix(c(1,1,0,0,0,0,
              0,0,1,1,0,0,
              0,0,0,0,1,1),6,3)
Z

Va <- 8
Ve <- 6

R <- diag(6)*Ve
R
Vo <- Va/4
G <- diag(3)*Vo
G

V <- Z %*% G %*% t(Z) + R
V

# BLUE
b <- solve(t(X)%*%solve(V)%*%X) %*% (t(X)%*%solve(V)%*%y)
b

# BLUP  
u <- G %*% t(Z) %*%solve(V) %*%(y-(X%*%b))
u

   # řešení pomocí soustavy normálních rovnic
XX <- t(X)%*%solve(R)%*%X
XX

XZ <- t(X)%*%solve(R)%*%Z
XZ

ZX <- t(Z)%*%solve(R)%*%X
ZX

ZZG <- t(Z)%*%solve(R)%*%Z + solve(G)
ZZG

LS1 <- cbind(XX,XZ)
LS2 <- cbind(ZX,ZZG)
LS <- rbind(LS1,LS2)
LS

Xy <- t(X)%*%solve(R)%*%y
Zy <- t(Z)%*%solve(R)%*%y
PS <- rbind(Xy,Zy)
PS

bu <- solve(LS) %*% PS
bu

### 2)  výpočet matice aditivní příbuznosti A pomocí balíčku pedigreemm) z příkladu rodokmen.pdf

library(pedigreemm)

ped <- pedigree(sire = c(NA,NA,1,1,4,5),
                dam = c(NA,NA,2,NA,3,2), label = 1:6)
ped

A <- getA(ped)
A

A <- matrix(A,6,6)
A



# 3)  animal model - výpočet plemenné hodnoty pro všechny jedince z databáze

ped <- pedigree(sire = c(NA,NA,NA,1,3,3),
                dam = c(NA,NA,NA,2,2,4), label=1:6)
ped

A <- getA(ped)
A <- matrix(A,6,6)
A

y <- matrix(c(200,170,180),3,1)

X <- matrix(c(1,1,1),3,1)
X

Za <- diag(0,3)
Za
Zb <- diag(1,3)
Zb
Z <- cbind(Za,Zb)
Z

h2 <- 1/3
k <- (1-h2)/h2
k

XX <- t(X)%*%X
XX

XZ <- t(X)%*%Z
XZ

ZX <- t(Z)%*%X
ZX

ZZAk <- t(Z)%*%Z + solve(A)*k
ZZAk

LS <- rbind((cbind(XX,XZ)),(cbind(ZX,ZZAk)))
LS

Xy <- t(X)%*%y
Xy
Zy <- t(Z)%*%y
Zy

PS <- rbind(Xy,Zy)
PS

bu <- solve(LS) %*% PS
bu