1# These functions are 2# Copyright (C) 1998-2021 T.W. Yee, University of Auckland. 3# All rights reserved. 4 5 6 7 8 psv2magic <- 9 function(x.VLM, constraints, spar.vlm, sm.osps.list) { 10 11 12 13 14 colperm <- function(x, from, to) { 15 16 ncx <- ncol(x) 17 if (length(from) != length(to) || 18 any(from != round(from)) || 19 any(from < 1 | ncx < from) || 20 any(duplicated(from)) || 21 any(sort(from) != sort(to))) 22 stop("invalid column permutation indices") 23 perm <- seq_len(ncx) 24 perm[to] <- perm[from] 25 x[, perm] 26 } 27 28 29 30 assignx <- sm.osps.list$assignx 31 nassignx <- names(assignx) 32 indexterms <- sm.osps.list$indexterms 33 which.X.sm.osps <- sm.osps.list$which.X.sm.osps 34 term.labels <- sm.osps.list$term.labels 35 ncol.X.sm.osps <- sapply(which.X.sm.osps, length) 36 ncolHlist.model <- unlist(lapply(constraints, ncol)) 37 38 39 ncolHlist.new <- ncolHlist.model 40 if (names(constraints)[[1]] == "(Intercept)") { 41 ncolHlist.new <- ncolHlist.new[-1] 42 nassignx <- nassignx[-1] 43 } 44 45 46 ncol.H.ps <- ncolHlist.new[indexterms] 47 num.osps.terms <- length(which.X.sm.osps) 48 49 50 ncol.allterms <- sapply(assignx, length) 51 52 ncol.model <- if (names(constraints)[[1]] == "(Intercept)") 53 ncol.allterms[-1] else ncol.allterms 54 jay <- 0 55 jjoffset <- if (names(constraints)[[1]] == "(Intercept)") 56 ncolHlist.model[1] else 0 57 perm.list <- list() 58 for (ii in seq_along(term.labels)) { 59 if (indexterms[ii]) { 60 jay <- jay + 1 61 perm.list[[jay]] <- 62 matrix(jjoffset + 1:(ncol.X.sm.osps[jay] * ncol.H.ps[jay]), 63 nrow = ncol.X.sm.osps[jay], # Redundant really 64 ncol = ncol.H.ps[jay], byrow = TRUE) 65 jjoffset <- jjoffset + ncol.H.ps[[jay]] * ncol.X.sm.osps[[jay]] 66 } else { 67 jjoffset <- jjoffset + ncolHlist.new[ii] * ncol.model[ii] 68 } 69 } # for ii 70 vindex.min <- sapply(perm.list, min) # function(x) min(x) 71 vindex.max <- sapply(perm.list, max) # function(x) max(x) 72 oo1 <- vector("list", length(ncol.H.ps)) # list() 73 for (ii in seq_along(ncol.H.ps)) { 74 oo1[[ii]] <- seq.int(vindex.min[ii], vindex.max[ii]) 75 } 76 ooo <- unlist(oo1, use.names = FALSE) # do.call("c", oo1) 77 ppp <- unlist(perm.list, use.names = FALSE) # do.call("c", perm.list) 78 79 80 OFF.list <- vector("list", num.osps.terms) # list() 81 for (ii in 1:num.osps.terms) { 82 index <- 0 83 OFF.list[[ii]] <- numeric() 84 for (jay in 1:(ncol.H.ps[ii])) { 85 OFF.list[[ii]][jay] <- vindex.min[ii] + index 86 index <- ncol.X.sm.osps[ii] * jay 87 } 88 } 89 90 91 list(x.VLM.new = if (identical(ppp, ooo)) x.VLM else 92 colperm(x.VLM, ppp, ooo), 93 sp = unlist(spar.vlm), 94 S.arg = rep(sm.osps.list$S.arg, ncol.H.ps), # Argument 'S' of magic() 95 OFF = unlist(OFF.list)) 96} # psv2magic 97 98 99