1
$\begingroup$

I am working in R, and am trying to generate values of

$$ logit^{-1}(\alpha X_1+\beta X_2) $$

with $\alpha,\beta$ such that $logit^{-1}(\alpha X_1+\beta X_2)$ has a mean of 0.4 with as small a variance as possible. I have that $X_1 \sim Bern(0.5)$ and that $X_2\sim N(0.5,0.01)$.

My code is:

X1 <- rbinom(100, size = 1, prob = 0.5) 
X2 <- rnorm(100, mean=0.5, sd=0.1)    
logit.inverse <- function(x) exp(x)/(1+exp(x))
alpha <- -5
beta  <- 5
mean(logit.inverse(alpha*X1+beta*X2))
var(logit.inverse(alpha*X1+beta*X2))

The results above are $0.5175315$ for the mean, and $0.1776111$ for the variance.

I am wondering if there is an optimization technique to make it so that the mean is $0.4$ or $0.3$ with as small variance as possible? Thanks!

$\endgroup$
1
$\begingroup$

The next best thing to a simple formula for an expectation of a random variable is a zillion random samples for that random variable. Below is some R code to estimate the values of $\alpha$ and $\beta$ that have expectation (approximately) 0.4 and (approximate) minimum variance.

First a large number of random samples are taken and a contour plot of the expectations associated with $\alpha$ and $\beta$ to get a lay of the land:

# Generate a large random sample
  nsim <- 10000
  x1 <- rbinom(nsim, size = 1, prob = 0.5) 
  x2 <- rnorm(nsim, mean=0.5, sd=0.1)    

# Function to estimate mean and variance of 
#  1 - 1/(1+exp(alpha*X1 + beta*X2))
  logit.inverse <- function(a, b) {
    data.frame(mean = mean(1 - 1/(exp(a*x1 + b*x2))),
      var = var(1 - 1/(exp(a*x1 + b*x2))))
  }

# Create a grid of estimated expectations 
  xmin <- -2
  xmax <- 2
  ymin <- -1
  ymax <- 3
  n <- 100 # Number of grid cells in each direction
  z <- matrix(rep(NA,(n+1)^2), nrow=n+1)
  for (i in 1:(n+1)) {
      a <- xmin + (xmax-xmin)*(i-1)/n
      for (j in 0:n) {
          b <- ymin + (ymax-ymin)*(j-1)/n
          z[i,j] <- logit.inverse(a,b)$mean
      }                                                
  }

# Create a contour plot to see what's what
  par(ask=TRUE)
  contour(xmin + (xmax-xmin)*c(0:n)/n, ymin + (ymax-ymin)*c(0:n)/n, z,
    levels = c(0,0.1,0.2,0.3,0.4,0.5,0.6,0.7), xlab="alpha", ylab="beta",
    labcex=1, las=1, main="E(1-1/(1+exp(alpha*X1+beta*X2)))")

# Get contour line for the expectation = 0.4
  k <- contourLines(xmin + (xmax-xmin)*c(0:n)/n, ymin + (ymax-ymin)*c(0:n)/n, z,
    levels = 0.4)[[1]]

# Determine alpha and beta associated with the minimum variance
  v <- rep(NA,length(k$x))
  for (i in 1:length(k$x)) {
      v[i] <- logit.inverse(k$x[i],k$y[i])$var 
  }

# Add in the "optimal" point
  points(k$x[v==min(v)], k$y[v==min(v)], pch=16, col="red", cex=2)

Contour plot of expectations

# Summarize results
  cat("alpha =", k$x[v==min(v)],"\n")
  cat("beta =", k$y[v==min(v)],"\n")
  cat("Minimum variance =", min(v),"\n")

alpha = 0.04 beta = 0.995488 Minimum variance = 0.00372371

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.