simulate_burninthining_multiple: simulate_burninthining_multiple
In isci1102/ERPM:

View source: R/functions_burninthining.R

1	simulate_burninthining_multiple(partitions, presence.tables, theta, nodes, effects, objects, num.steps, neighborhood, sizes.allowed, sizes.simulated, max.thining)

`partitions`	Observed partitions.
`presence.tables`	Matrix indicating which actors were present for each observations (mandatory).
`theta`	Initial model parameters
`nodes`	Data frame containing the nodes.
`effects`	Effects or sufficient statistics. A list with a vector "names" and a vector "objects".
`objects`	Objects used for statistics calculation. A list with a vector "name" and a vector "object".
`num.steps`	Number of samples wanted.
`neighborhood`	Way of choosing partitions: probability vector (probability actors swap, probability merge/division, probability single actor move).
`sizes.allowed`	Vector of group sizes allowed in sampling.
`sizes.simulated`	Vector of group sizes allowed in the Markov chain but not necessraily sampled.
`max.thining`	Where to stop adding thining.

##---- Should be DIRECTLY executable !! ----
##-- ==>  Define data, use random,
##--	or do  help(data=index)  for the standard data sets.

## The function is currently defined as
function (partitions, presence.tables, theta, nodes, effects,
    objects, num.steps, neighborhood, sizes.allowed, sizes.simulated,
    max.thining)
{
    num.effects <- length(effects$names)
    chain <- draw_Metropolis_multiple(theta, partitions, presence.tables,
        nodes, effects, objects, 1, 1, num.steps, neighborhood,
        sizes.allowed, sizes.simulated, return.all.partitions = F)
    ok <- F
    current.draws <- chain$draws
    current.last <- chain$last.partitions
    current.thining <- 1
    allautocors <- c()
    while (!ok) {
        draws <- current.draws[seq(1, nrow(current.draws), current.thining),
            ]
        autocors <- rep(0, num.effects)
        for (e in 1:num.effects) {
            autocors[e] <- cor(draws[1:(num.steps - 1), e], draws[2:num.steps,
                e])
        }
        if (current.thining%%50 == 0) {
            print("thining")
            print(current.thining)
            print("autocorrelations")
            print(autocors)
        }
        allautocors <- rbind(allautocors, autocors)
        if (current.thining >= max.thining) {
            ok <- T
        }
        else {
            new.chain <- draw_Metropolis_multiple(theta, current.last,
                presence.tables, nodes, effects, objects, 1,
                1, num.steps, neighborhood, sizes.allowed, sizes.simulated,
                return.all.partitions = F)
            current.draws <- rbind(current.draws, new.chain$draws)
            current.last <- new.chain$last.partitions
            current.thining <- current.thining + 1
        }
    }
    allmeans <- matrix(0, nrow(draws), num.effects)
    allmeans[1, ] <- draws[1, ]
    for (d in 2:nrow(draws)) {
        allmeans[d, ] <- colMeans(draws[1:d, ])
    }
    smoothedautocors <- allautocors
    smoothedmeans <- allmeans
    for (eff in 1:num.effects) {
        lo <- loess(y ~ x, data.frame(x = 1:max.thining, y = allautocors[,
            eff]))
        smoothedautocors[!is.na(smoothedautocors[, eff]), eff] <- lo$fitted
        lo <- loess(y ~ x, data.frame(x = 1:nrow(draws), y = allmeans[,
            eff]))
        smoothedmeans[!is.na(smoothedmeans[, eff]), eff] <- lo$fitted
    }
    return(list(draws = current.draws, autocorrelations = allautocors,
        moving.means = allmeans, autocorrelations.smoothed = smoothedautocors,
        moving.means.smoothed = smoothedmeans))
  }