Skip to content

Instantly share code, notes, and snippets.

@tjmahr
Last active April 17, 2024 16:51
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save tjmahr/1a0298457dfb9a4595b7730e01edb59f to your computer and use it in GitHub Desktop.
Save tjmahr/1a0298457dfb9a4595b7730e01edb59f to your computer and use it in GitHub Desktop.
x <- rnorm(250, 5, 2)
i <- sample.int(250, 20)

ranks <- rank(x)[i]
percentiles <- ppoints(250)[ranks]
points <- x[i]

f <- function(args, points, percentiles) {
  test_quantiles <- qnorm(percentiles, args[1], args[2])
  sum((test_quantiles - points) ^ 2)
}
l <- optim(c(1, 1), f, points = points, percentiles = percentiles)
l
#> $par
#> [1] 5.023273 1.861110
#> 
#> $value
#> [1] 0.2986559
#> 
#> $counts
#> function gradient 
#>       61       NA 
#> 
#> $convergence
#> [1] 0
#> 
#> $message
#> NULL


library(ggplot2)
theoretical <- tibble::tibble(
  x = ppoints(250),
  y = qnorm(x, l$par[1], l$par[2])
)
subsample <- tibble::tibble(
  x = percentiles,
  y = points
)
ggplot(theoretical) + 
  aes(x = x, y = y) +
  geom_line() +
  geom_point(data = subsample, size = 3, color = "orange")

x <- rexp(250, 1 / 2500)
i <- sample.int(250, 20)

ranks <- rank(x)[i]
percentiles <- ppoints(250)[ranks]
points <- x[i]

fexp <- function(args, points, percentiles) {
  test_quantiles <- qexp(percentiles, args[1])
  sum((test_quantiles - points) ^ 2)
}

l <- optim(1, fexp, points = points, percentiles = percentiles)
#> Warning in optim(1, fexp, points = points, percentiles = percentiles): one-dimensional optimization by Nelder-Mead is unreliable:
#> use "Brent" or optimize() directly
#> Warning in qexp(percentiles, args[1]): NaNs produced
#> Warning in qexp(percentiles, args[1]): NaNs produced
#> Warning in qexp(percentiles, args[1]): NaNs produced
#> Warning in qexp(percentiles, args[1]): NaNs produced
#> Warning in qexp(percentiles, args[1]): NaNs produced
#> Warning in qexp(percentiles, args[1]): NaNs produced
#> Warning in qexp(percentiles, args[1]): NaNs produced
#> Warning in qexp(percentiles, args[1]): NaNs produced
#> Warning in qexp(percentiles, args[1]): NaNs produced
#> Warning in qexp(percentiles, args[1]): NaNs produced
#> Warning in qexp(percentiles, args[1]): NaNs produced
#> Warning in qexp(percentiles, args[1]): NaNs produced
#> Warning in qexp(percentiles, args[1]): NaNs produced
#> Warning in qexp(percentiles, args[1]): NaNs produced
#> Warning in qexp(percentiles, args[1]): NaNs produced
#> Warning in qexp(percentiles, args[1]): NaNs produced
#> Warning in qexp(percentiles, args[1]): NaNs produced
#> Warning in qexp(percentiles, args[1]): NaNs produced
#> Warning in qexp(percentiles, args[1]): NaNs produced
#> Warning in qexp(percentiles, args[1]): NaNs produced
#> Warning in qexp(percentiles, args[1]): NaNs produced
l
#> $par
#> [1] 0.0004137039
#> 
#> $value
#> [1] 838223.3
#> 
#> $counts
#> function gradient 
#>       70       NA 
#> 
#> $convergence
#> [1] 0
#> 
#> $message
#> NULL

library(ggplot2)
theoretical <- tibble::tibble(
  x = ppoints(250),
  y = qexp(x, l$par[1])
)
subsample <- tibble::tibble(
  x = percentiles,
  y = points
)
ggplot(theoretical) + 
  aes(x = x, y = y) +
  geom_line() +
  geom_point(data = subsample, size = 3, color = "orange")

Created on 2024-04-17 with reprex v2.1.0

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