和加权瑞利分布

Question

我根据帕累托分布函数使用nrow=1000（个人）和ncol=100（天）模拟数据的步长：

set.seed(10)
sim_data <- replicate(100, VGAM::rpareto(1000, shape=7, scale=500))
sim_data <- as.data.frame(sim_data)
set.seed(10)
sim_data[,1:50] <- sim_data[50]*(-1) ##assign directionality
sim_data_directions <- as.data.frame(sim_data)
##randomize columns
sim_data <- sim_data_directions[,sample(ncol(sim_data_directions))] 
x <- c(1:31) ## for variable 'days' to associate to each step length
set.seed(10)
df_sim <- cbind(sim_data, t(apply(sim_data,1, function(x) {
  i1 <- sample(seq_along(x), 1)
  out <- sum(sample(x, i1))
  c(days = i1, step_lengths = out)}
))) ## create step lengths
##adding weights for each level in variable days
df_sim <- as.data.frame(dplyr::add_count(df_sim, days))

使用此数据集df_sim，其中包含步长的模拟值、与每个步长相关的时间（以天为单位）和权重（每个时间变量的值数（以天为单位），我想总结分布，使用一个瑞利分布函数，其中每个级别的天数的分布是加权的，如下所示：

rayleigh_distr <- sum(n*function (x) x*exp(-1*(x/2*sigma)^2)/sigma^2)

n 是权重。 如何根据权重总结每天的分布？

Answer 1

首先，您的 Rayleigh PDF 似乎有错误。应该是：

x*exp(-(x/sigma)^2/2)/sigma^2

您似乎希望将 rayleigh_distr 作为一个函数，返回带有权重 n 和比例参数 sigma 的 Rayleigh mixture distribution 的 PDF。如果是这样的话，那就是（对于标量x）：

rayleigh_distr <- function(x, sigma, n) sum(n*x*exp(-(x/sigma)^2/2)/sigma^2)

其中sigma 和n 是长度相等的向量。

如果您想为x 传递一个向量（不一定与sigma 和n 的长度相同），那么这应该可以：

library(Rfast)
rayleigh_distr <- function(x, sigma, n) colsums(n*eachrow(exp(-(outer(1/sigma, x))^2/2), x)/sigma^2)

它将返回一个与x 长度相同的向量。

更新：

CDF 将是：

rayleigh_cdf <- function(x, sigma, n) colSums(n*(1 - exp(-(outer(1/sigma, x))^2/2)))