Two-dimensional Monte Carlo: Difference between revisions

From Opasnet
Jump to navigation Jump to search
(→‎Rationale: updated to enable additional indices using info ovariable)
(→‎Rationale: strength added to parameters)
 
(3 intermediate revisions by the same user not shown)
Line 12: Line 12:
*  N2 = 1000, # Number of iterations in the new Iter
*  N2 = 1000, # Number of iterations in the new Iter
*  run2d = TRUE, # Should the mc2d function be used or not?
*  run2d = TRUE, # Should the mc2d function be used or not?
*  newmarginals = c("Gender", "Ages", "Country"), # Names of columns that are non-marginals but will be sampled enough to become marginals {{attack|# |The function will produce an ovariable that correctly has these indices as marginals. However, if the function is used within an ovariable formula (which is typically the case), the marginal status is in the end inherited from parents and they are re-converted to non-marginals. This should be fixed somehow.|--[[User:Jouni|Jouni]] ([[User talk:Jouni|talk]]) 15:11, 11 June 2017 (UTC)}}
*  info = 1, # An ovariable that may add new indices to the ovariable to be converted. If none, use 1.
*  newmarginals = c("Gender", "Ages", "Country"), # Names of columns that are non-marginals but will be sampled enough to become marginals. The function will produce an ovariable that correctly has these indices as marginals. However, if the function is used within an ovariable formula (which is typically the case), the marginal status is in the end inherited from parents and they may be re-converted to non-marginals. If this happens, the marginal status has to be updated in the assessment model code on case by case basis. Any automatic solution would violate the inheritance rules.
*  method = "bootstrap", # which method to use for 2D Monte Carlo? Currently bootsrap is the only option.
*  method = "bootstrap", # which method to use for 2D Monte Carlo? Currently bootsrap is the only option.
*  fun = mean # Function for aggregating the first Iter dimension.
*  fun = mean # Function for aggregating the first Iter dimension.
Line 20: Line 21:
  objects.latest("Op_en7805", code_name = "mc2d")
  objects.latest("Op_en7805", code_name = "mc2d")


== Rationale ==


* 22.11.2017 Version where run2d==FALSE does nothing. [http://en.opasnet.org/en-opwiki/index.php?title=Special:RTools&id=QUs11wbaRAVtNLgs]


== Rationale ==
<rcode name="mc2d" label="Initiate function mc2d" embed=1>
 
<rcode name="mc2d" embed=1>
#This is code Op_en7805/mc2d on page [[Two-dimensional Monte Carlo]]
#This is code Op_en7805/mc2d on page [[Two-dimensional Monte Carlo]]
library(OpasnetUtils)
library(OpasnetUtils)
Line 34: Line 35:
### Sitten näistä arvotaan bootstrapillä uudet iteraatiot jotka aggregoidaan populaatiotasolle.
### Sitten näistä arvotaan bootstrapillä uudet iteraatiot jotka aggregoidaan populaatiotasolle.
### Tässä versiossa on aina sama RR-Iter-kombinaatio.
### Tässä versiossa on aina sama RR-Iter-kombinaatio.
### strength: Estimated number of actual observations (number of simulations does not reflect this)
## Ensin muutetaan yksilötason Iter Iter2:ksi, jolloin uusi indeksi paisuttaa ovariablen N-kertaiseksi.
## Ensin muutetaan yksilötason Iter Iter2:ksi, jolloin uusi indeksi paisuttaa ovariablen N-kertaiseksi.
### Sitten lasketaan ja lopuksi aggregoidaan populaatiotasolle.
### Sitten lasketaan ja lopuksi aggregoidaan populaatiotasolle.
Line 40: Line 42:
# Boostrap-versio:
# Boostrap-versio:


# Paramter list. Note: this is not stored, you have to define it in the model code.
mc2d <- function(ova, mc2dpar = NULL) {
 
   if(is.null(mc2dpar)) if(exists("mc2dparam")) mc2dpar <- mc2dparam else stop("Parameter list mc2dparam missing!\n")
mc2dparam<- list(
   if(mc2dpar$run2d) {
  N2 = 1000, # Number of iterations in the new Iter
    ova <- ova * mc2dpar$info
  run2d = TRUE, # Should the mc2d function be used or not?
    require(reshape2)
   info = 1, # Ovariable that contains additional indices, e.g. newmarginals. If none, use 1.
    marg <- setdiff(c(colnames(ova@output)[ova@marginal], mc2dpar$newmarginals), "Iter")
  newmarginals = c("Gender", "Ages", "Country"), # Names of columns that are non-marginals but should be sampled enough to become marginals
    out <- aggregate(
  method = "bootstrap", # which method to use for 2D Monte Carlo? Currently bootsrap is the only option.
      result(ova),
  fun = mean # Function for aggregating the first Iter dimension.
      by = ova@output[colnames(ova@output) %in% marg],
)
      FUN = function(x) {  
 
        strength <- if(is.null(mc2dpar$strength)) length(x) else mc2dpar$strength
mc2d <- function(ova) {
        apply(
  if(!exists("mc2dparam")) stop("Parameter list mc2dparam missing!\n")
          array(
   ova <- ova * mc2dparam$info
            as.numeric(sample(as.character(x), strength * mc2dpar$N2, replace=TRUE)),
  require(reshape2)
            #Numeric conversion is needed to prevent x from being interpreted as number of choices.
  marg <- setdiff(c(colnames(ova@output)[ova@marginal], mc2dparam$newmarginals), "Iter")
            dim = c(strength,mc2dpar$N)
  out <- aggregate(
          ),
    result(ova),
          MARGIN=2, FUN = mc2dpar$fun
    by = ova@output[colnames(ova@output) %in% marg],
        )
    FUN = function(x) {  
      }
      apply(
    )
        array(
    temp <- melt(out[[length(out)]])
          as.numeric(sample(as.character(x), length(x)*mc2dparam$N2, replace=TRUE)),
    out[[length(out)]] <- 1:nrow(out)
          #Numeric conversion is needed to prevent x from being interpreted as number of choices.
    colnames(temp) <- c("Nrow","Iter","Result")
          dim = c(length(x),mc2dparam$N)
    out <- merge(out, temp, by.x = "x", by.y="Nrow")
        ),
    out$x <- NULL
        MARGIN=2, FUN = mc2dparam$fun
    out <- Ovariable(
      )
      output = out,  
    }
      marginal = colnames(out) %in% c(marg, "Iter")
  )
    )
  temp <- melt(out[[length(out)]])
  } else {
  out[[length(out)]] <- 1:nrow(out)
    out <- ova
  colnames(temp) <- c("Nrow","Iter","Result")
   }
  out <- merge(out, temp, by.x = "x", by.y="Nrow")
  out$x <- NULL
  out <- Ovariable(
    output = out,  
    marginal = colnames(out) %in% c(marg, "Iter")
   )
   return(out)
   return(out)
}
}
Line 84: Line 80:
objects.store(mc2d)
objects.store(mc2d)
cat("Function mc2d stored.\n")
cat("Function mc2d stored.\n")
</rcode>
<rcode name="mc2dparam" label="Initiate mc2dparam" embed=1>
#This is code Op_en7805/mc2dparam on page [[Two-dimensional Monte Carlo]]
library(OpasnetUtils)
# Parameter list.
mc2dparam<- list(
#  N2 = 1000, # Number of iterations in the new Iter
  run2d = FALSE # Should the mc2d function be used or not?
#  info = 1, # Ovariable that contains additional indices, e.g. newmarginals. If none, use 1.
#  newmarginals = c("Gender", "Ages", "Country"), # Names of columns that are non-marginals but should be sampled enough to become marginals
#  method = "bootstrap", # which method to use for 2D Monte Carlo? Currently bootsrap is the only option.
#  fun = mean # Function for aggregating the first Iter dimension.
)
objects.store(mc2dparam)
cat("Function mc2dparam stored.\n")
</rcode>
</rcode>


Latest revision as of 13:35, 17 June 2019



Question

How to perform two-dimensional Monte Carlo in Opasnet?

Answer

Use function mc2d to perform two-dimensional Monte Carlo. The function samples the current ovariable results by bootstrapping, applies an aggregate function to the samples, and then produces a new Iter index location for each sample. The function requires a parameter list mc2dparam, which contains the following parameters (with some example values):

  • N2 = 1000, # Number of iterations in the new Iter
  • run2d = TRUE, # Should the mc2d function be used or not?
  • info = 1, # An ovariable that may add new indices to the ovariable to be converted. If none, use 1.
  • newmarginals = c("Gender", "Ages", "Country"), # Names of columns that are non-marginals but will be sampled enough to become marginals. The function will produce an ovariable that correctly has these indices as marginals. However, if the function is used within an ovariable formula (which is typically the case), the marginal status is in the end inherited from parents and they may be re-converted to non-marginals. If this happens, the marginal status has to be updated in the assessment model code on case by case basis. Any automatic solution would violate the inheritance rules.
  • method = "bootstrap", # which method to use for 2D Monte Carlo? Currently bootsrap is the only option.
  • fun = mean # Function for aggregating the first Iter dimension.

You can call the function by using code 

objects.latest("Op_en7805", code_name = "mc2d")

Rationale

  • 22.11.2017 Version where run2d==FALSE does nothing. [1]

+ Show code

+ Show code

See also

Retrieved from "https://en.opasnet.org/index.php?title=Two-dimensional_Monte_Carlo&oldid=42864"