# R/Matrixmult.R In RandVar: Implementation of Random Variables

```## "Matrixmultiplication"
setMethod("%*%", signature(x = "matrix", y = "EuclRandVariable"),
function(x, y){
if(!is.numeric(x))
stop("x is no numeric matrix")
if(ncol(x) != dimension(y))
stop("the number of columns of x != dimension of y")

fct1 <- NULL;  e <- NULL
dimn <- y@Range@dimension
map <- vector("list", nrow(x))
for(i in 1:nrow(x)){
if(length(y) == 1){
if(!identical(all.equal(x[i,1:dimn], numeric(dimn),
tolerance = .Machine\$double.eps^0.75), TRUE)){
map[[i]] <- function(x){ f1 <- fct1; e %*% f1(x) }
body(map[[i]]) <- substitute({ f1 <- fct1; e %*% f1(x) },
list(fct1 = y@Map[[1]], e = x[i,1:dimn]))
}else{
map[[i]] <- function(x){ numeric(1) }
}
}else{
if(!identical(all.equal(x[i, 1:dimn], numeric(dimn),
tolerance = .Machine\$double.eps^0.75), TRUE)){
fct <- function(x){ f1 <- fct1; e %*% f1(x) }
body(fct) <- substitute({ f1 <- fct1; e %*% f1(x) },
list(fct1 = y@Map[[1]], e = x[i,1:dimn]))
}else{
fct <- function(x) { numeric(1) }
}
for(j in 2:length(y)){
if(!identical(all.equal(x[i, ((j-1)*dimn + 1):(j*dimn)], numeric(dimn),
tolerance = .Machine\$double.eps^0.75), TRUE)){
body(fct) <- substitute({ f1 <- fct1; f2 <- fct2; f1(x) + e %*% f2(x)},
list(fct1 = fct, fct2 = y@Map[[j]],
e = x[i,((j-1)*dimn + 1):(j*dimn)]))
}
}
map[[i]] <- fct
}
}

R <- new("EuclRandMatrix")
R@Map <- map
R@Dim <- as.integer(c(nrow(x), 1))
R@Domain <- y@Domain
R@Range <- Reals()

return(R)
})
setMethod("%*%", signature(x = "matrix", y = "EuclRandMatrix"),
function(x, y){
dx <- dim(x)
dy <- y@Dim
dimn <- y@Range@dimension
if(dx[2] != dy[1]*dimn)
stop("non-conformable arguments")
map <- vector("list", dx[1]*dy[2])

for(i in 1:dy[2])
for(j in 1:dx[1])
map[[(i-1)*dx[1] + j]] <- Map(t(x[j,]) %*% y[,i])[[1]]

y@Map <- map
y@Dim <- c(dx[1], dy[2])
y@Range <- Reals()

return(y)
})
setMethod("%*%", signature(x = "matrix", y = "EuclRandVarList"),
function(x, y){
if(!is.numeric(x))
stop("x is no numeric matrix")
if(ncol(x) != dimension(y))
stop("the number of columns of x != dimension of y")

fct1 <- NULL;  e <- NULL
map <- vector("list", nrow(x))
for(i in 1:nrow(x)){
if(numberOfMaps(y) == 1){
dimn <- y[[1]]@Range@dimension
if(!identical(all.equal(x[i,1:dimn], numeric(dimn),
tolerance = .Machine\$double.eps^0.75), TRUE)){
map[[i]] <- function(x){ f1 <- fct1; e %*% f1(x) }
body(map[[i]]) <- substitute({ f1 <- fct1; e %*% f1(x) },
list(fct1 = y[[1]]@Map[[1]], e = x[i,1:dimn]))
}else{
map[[i]] <- function(x){ numeric(1) }
}
}else{
comp <- 0
ind2 <- 0
for(k in 1:length(y)){
dimn <- y[[k]]@Range@dimension
for(j in 1:length(y[[k]])){
comp <- comp + 1
ind1 <- ind2 + 1
ind2 <- ind1 + dimn - 1
if(comp == 1){
if(!identical(all.equal(x[i, 1:dimn], numeric(dimn),
tolerance = .Machine\$double.eps^0.75), TRUE)){
fct <- function(x){ f1 <- fct1; e %*% f1(x) }
body(fct) <- substitute({ f1 <- fct1; e %*% f1(x) },
list(fct1 = y[[k]]@Map[[1]], e = x[i,1:dimn]))
}else{
fct <- function(x) { numeric(1) }
}
next
}
if(!identical(all.equal(x[i, ind1:ind2], numeric(dimn),
tolerance = .Machine\$double.eps^0.75), TRUE)){
body(fct) <- substitute({ f1 <- fct1; f2 <- fct2; f1(x) + e %*% f2(x)},
list(fct1 = fct, fct2 = y[[k]]@Map[[j]],
e = x[i, ind1:ind2]))
}
}
}
map[[i]] <- fct
}
}

R <- new("EuclRandMatrix")
R@Map <- map
R@Dim <- as.integer(c(nrow(x), 1))
R@Domain <- y[[1]]@Domain
R@Range <- Reals()

return(R)
})
setMethod("%*%", signature(x = "numeric", y = "EuclRandMatrix"),
function(x, y){ t(x) %*% y })

setMethod("%*%", signature(x = "numeric", y = "EuclRandVariable"),
function(x, y){ t(x) %*% as(y, "EuclRandMatrix") })

setMethod("%*%", signature(x = "EuclRandVariable", y = "matrix"),
function(x, y){ t(y %*% x) })

setMethod("%*%", signature(x = "EuclRandMatrix", y = "matrix"),
function(x, y){ t(t(y) %*% t(x)) })

setMethod("%*%", signature(x = "EuclRandVarList", y = "matrix"),
function(x, y){ t(y %*% x) })

setMethod("%*%", signature(x = "EuclRandMatrix", y = "numeric"),
function(x, y){ x %*% as.matrix(y) })

setMethod("%*%", signature(x = "EuclRandVariable", y = "numeric"),
function(x, y){ t(x) %*% as.matrix(y) })

setMethod("%*%", signature(x = "EuclRandVariable", y = "EuclRandVariable"),
function(x, y){
nrvalues <- length(x)
if(nrvalues != length(y))
stop("non-conformable arguments")
if(!compatibleDomains(x, y))
stop("the domains of the two random variables are not compatible")

fct1 <- NULL;  fct2 <- NULL
fct <- function(x){ f1 <- fct1; f2 <- fct2; f1(x) %*% f2(x) }
body(fct) <- substitute({ f1 <- fct1; f2 <- fct2; f1(x) %*% f2(x) },
list(fct1 = x@Map[[1]], fct2 = y@Map[[1]]))

if(nrvalues > 1)
for(i in 2:nrvalues)
body(fct) <- substitute({ f <- fct; f1 <- fct1; f2 <- fct2; f(x) + f1(x) %*% f2(x) },
list(fct = fct, fct1 = x@Map[[i]], fct2 = y@Map[[i]]))

if(identical(x@Domain, y@Domain))
Dom <- x@Domain
else
if(is(x@Domain, class(y@Domain)[1]))
Dom <- x@Domain
else
Dom <- y@Domain

R <- new("EuclRandMatrix")
R@Map <- list(fct)
R@Dim <- as.integer(c(1, 1))
R@Domain <- Dom
R@Range <- Reals()

return(R)
})
setMethod("%*%", signature(x = "EuclRandMatrix", y = "EuclRandMatrix"),
function(x, y){
dx <- x@Dim
dy <- y@Dim
if(dx[2] != dy[1])
stop("non-conformable arguments")
if(!compatibleDomains(x, y))
stop("the domains of the two random variables are not compatible")
if(x@Range@dimension != y@Range@dimension)
stop("the two random matrices have different ranges")
map <- vector("list", dx[1]*dy[2])

for(i in 1:dy[2])
for(j in 1:dx[1])
map[[(i-1)*dx[1] + j]] <- Map(x[j,] %*% y[,i])[[1]]

x@Map <- map
x@Dim <- c(dx[1], dy[2])
x@Range <- EuclideanSpace(dimension = x@Range@dimension*y@Range@dimension)

return(x)
})
setMethod("%*%", signature(x = "EuclRandVariable", y = "EuclRandMatrix"),
function(x, y){
if(y@Dim[1] != 1)
stop("non-conformable arguments")
if(!compatibleDomains(x, y))
stop("the domains of the two random variables are not compatible")

R <- new("EuclRandMatrix")
R@Map <- x@Map
R@Dim <- as.integer(c(length(x), 1))
R@Domain <- x@Domain
R@Range <- x@Range

return(R %*% y)
})
setMethod("%*%", signature(x = "EuclRandMatrix", y = "EuclRandVariable"),
function(x, y){
nrvalues <- length(y)
if(nrvalues != x@Dim[2])
stop("non-conformable arguments")
if(!compatibleDomains(x, y))
stop("the domains of the two random variables are not compatible")

R <- new("EuclRandMatrix")
R@Map <- y@Map
R@Dim <- as.integer(c(nrvalues, 1))
R@Domain <- y@Domain
R@Range <- y@Range

return(x %*% R)
})
```

## Try the RandVar package in your browser

Any scripts or data that you put into this service are public.

RandVar documentation built on April 6, 2019, 3 a.m.