Skip to content

Instantly share code, notes, and snippets.

@thomvolker
Last active September 27, 2024 10:50
Show Gist options
  • Save thomvolker/a93eedd5ee9983620b00eecc57352593 to your computer and use it in GitHub Desktop.
Save thomvolker/a93eedd5ee9983620b00eecc57352593 to your computer and use it in GitHub Desktop.
Density ratio estimation is invariant to one-to-one and onto transformations

Density ratio estimation is invariant to one-to-one and onto transformations

Thom Benjamin Volker

At least in the univariate case, the multivariate case must be considered in a future gist still.

N <- 10000000
x1 <- rexp(N, 1)

hist(x1, freq = FALSE, breaks = 30)
curve(dexp(x, 1), add = TRUE, col = "red")

x2 <- rexp(N, 2)
hist(x2, freq = FALSE, breaks = 30)
curve(dexp(x, 2), add = TRUE, col = "blue")

x1log <- log(x1)
hist(x1log, freq = FALSE, breaks = 30)
curve(dexp(exp(x), 1) * exp(x), add = TRUE, col = "red")

x2log <- log(x2)
hist(x2log, freq = FALSE, breaks = 30)
curve(dexp(exp(x), 2) * exp(x), add = TRUE, col = "blue")

y <- 1
dexp(y, 1) / dexp(y, 2)
#> [1] 1.359141
dexp(exp(log(1)), 1) * exp(log(1)) / dexp(exp(log(1)), 2) * exp(log(1))
#> [1] 1.359141
dexp(2, 1) / dexp(2, 2)
#> [1] 3.694528
dexp(exp(log(2)), 1) * exp(log(2)) / (dexp(exp(log(2)), 2) * exp(log(2)))
#> [1] 3.694528


z1 <- rlnorm(N)
hist(z1, freq = FALSE, breaks = 30)
curve(dlnorm(x), add = TRUE, col = "red")

z2 <- rgamma(N, 2, 1)
hist(z2, freq = FALSE, breaks = 30)
curve(dgamma(x, 2, 1), add = TRUE, col = "blue")

y <- 1
dlnorm(1) / dgamma(1, 2, 1)
#> [1] 1.084438
dlnorm(exp(log(1))) * exp(log(1)) / (dgamma(exp(log(1)), 2, 1) * exp(log(1)))
#> [1] 1.084438

z1log <- log(z1)
hist(z1log, freq = FALSE, breaks = 30)
curve(dlnorm(exp(x)) * exp(x), add = TRUE, col = "red")

z2log <- log(z2)
hist(z2log, freq = FALSE, breaks = 30)
curve(dgamma(exp(x), 2, 1) * exp(x), add = TRUE, col = "blue")

y1 <- rnorm(2*N)
hist(y1, freq = FALSE, breaks = 30)
curve(dnorm(x), add = TRUE, col = "red")

y2 <- rt(2*N, 10)
hist(y2, freq = FALSE, breaks = 50)
curve(dt(x, 10), add = TRUE, col = "blue")

y1cb <- y1^3
hist(y1cb[y1cb > -5 & y1cb < 5], freq = FALSE, breaks = 100)
f1 <- \(y) 1 / (3 * sqrt(2 * pi)) * exp(-abs(y)^{2/3}/2) * abs(y)^{-2/3}
curve(f1(x), add = TRUE, col = "red", n=1000)

y2cb <- y2^3
hist(y2cb[y2cb > -40 & y2cb < 40], freq = FALSE, breaks = 200)
f2 <- \(y) gamma(5.5) / (3 * sqrt(10*pi) * gamma(5)) * abs(y)^{-2/3} *
  (1 + abs(y)^{2/3}/10)^{-11/2}
curve(f2(x), add = TRUE, col = "blue", n=1000)

dnorm(2) / dt(2, 10)
#> [1] 0.8829878
f1(2^3) / f2(2^3)
#> [1] 0.8829878

Created on 2024-09-27 with reprex v2.0.2

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment