Search code examples
rrandomestimation

random samples and distributions in R


I'm trying to build an estimator to compare the asymptotics of the GLS and OLS estimators. My idea is to try and see what happens at large samples, and with many of them.

Ideally, I would like to create a loop that would generate 6000 different random samples, of sizes 50 and 100 each, for different parameter values.

N=1000
n=c(50, 100)

#parameters
alpha0=1
beta0=1
gamma0=c(0, 0.1, 0.5)

alpha1=matrix(NA,N,6)
beta1=matrix(NA,N,6)
alpha2=matrix(NA,N,6)
beta2=matrix(NA,N,6)
alphaOLS=matrix(NA,N,6)
betaOLS=matrix(NA,N,6)

the different samples come from the combinations of gamma0 and n, which would equal 6 (times N) to get 6000. My first idea was to build a loop for the generation of the random samples the model I'm trying to work with is the following y_i=alpha+beta*x_i+u_i

u_i=e_i*h(x_i)^(1/2) and h(x)=exp(gamma0)

u <- list()


for (i in n) {
  for (k in gamma0) {
    x=rnorm(i,0,1)
    h=exp(gamma0[k]*x)
    e=rnorm(i,0,1)
    u[[i]] <- e*h^(1/2)
    
  }
  
}

The issue with this loop is that I'm only getting one random sample in x and e, and h is coming out as an empty matrix, and hence, u is also coming out empty. h here should be a matrix where the columns correspond to x* the different values of gamma0. e is supposed to be N(0,1) and u is meant to be the residual of the model

My ideal output should be get this loop to work, because from there on, I can sort my way around building an OLS and GLS estimator manually.

Thanks a lot!


Solution

  • This should work, here we use directly i and k instead of n and gamma0 in the loop.

    ### parameters
    N <- 1000
    n <- c(50, 100)
    
    alpha0 <- 1
    beta0 <- 1
    gamma0 <- c(0, 0.1, 0.5)
    
    ### Initiating objects for the loop
    u <- list()
    num_iter <- 0
    
    ### Looping
    for (i in n) {
      for (k in gamma0) {
        
        num_iter <- num_iter + 1
        
        x <- rnorm(i, 0, 1)
        h <- exp(k * x)
        
        e <- rnorm(i, 0, 1)
        u[[num_iter]] <- e * h^(1/2)
        names(u)[num_iter] <- paste("n:", i, ", gamma:", k, sep="", collapse=" ")
        
      }
    }
    
    ### Display results
    u
    

    Edit (based on your latest request)

    ### Parameters for each scenario
    N <- 1000
    n <- c(50, 100)
    gamma0 <- c(0, 0.1, 0.5)
    
    ### Initiating dataframes for each scenario
    for (i in n) {
      for (k in gamma0) {
    assign(paste("df_n", i, "gamma", k, sep=""), list())
      }
    }
    
    ### Defining a dataframe with n rows (sample size) and N columns(simulation) for each of the 6 scenarios
    df_n100gamma0 <- sapply(1:1000, function(x) rnorm(100, 0, 1) * exp(0 * rnorm(100, 0, 1))^(1/2))
    df_n100gamma0.1 <- sapply(1:1000, function(x) rnorm(100, 0, 1) * exp(0.1 * rnorm(100, 0, 1))^(1/2))
    df_n100gamma0.5 <- sapply(1:1000, function(x) rnorm(100, 0, 1) * exp(0.5 * rnorm(100, 0, 1))^(1/2))
    df_n50gamma0 <- sapply(1:1000, function(x) rnorm(50, 0, 1) * exp(0 * rnorm(50, 0, 1))^(1/2))
    df_n50gamma0.1 <- sapply(1:1000, function(x) rnorm(50, 0, 1) * exp(0.1 * rnorm(50, 0, 1))^(1/2))
    df_n50gamma0.5 <- sapply(1:1000, function(x) rnorm(50, 0, 1) * exp(0.5 * rnorm(50, 0, 1))^(1/2))
    
    ### Dimension of 1 dataframe
    dim(df_n100gamma0)
    
    ### Results
    # First simulation
    df_n100gamma0[, 1]
    
    # 57th simulation
    df_n100gamma0[, 57]