R中的对数转换（受其他约束）同时处理零

Question

我有一种情况，我需要对数据进行对数转换以使用它，但是矩阵中有零。除零外，我的矩阵也从Dirichlet分布中提取，这意味着该矩阵具有以下约束：所有列和必须相加为1。这是数据：

> q[1:10, 1:5]
            V1          V2          V3           V4           V5
1  0.534410243 0.009358740 0.011295181 0.2141751740 0.0030129254
2  0.026653603 0.372426720 0.447847534 0.0179177507 0.4072904477
3  0.193317915 0.003605024 0.003186611 0.4832114736 0.0007095471
4  0.111881585 0.000000000 0.000000000 0.2296213741 0.0119233461
5  0.089696570 0.591163629 0.509774416 0.0032542030 0.5535847030
6  0.007543558 0.000000000 0.000000000 0.0364907757 0.0013148362
7  0.004862942 0.000000000 0.002123909 0.0146682272 0.0004053690
8  0.009276195 0.011710457 0.014367894 0.0000000000 0.0000000000
9  0.006903171 0.004314528 0.011404455 0.0000000000 0.0126889937
10 0.015454219 0.007420903 0.000000000 0.0006610215 0.0090698319

注意q的所有列加起来为一

> colSums(q)[1:5]
V1 V2 V3 V4 V5 
 1  1  1  1  1 

我需要这样取log（q）：

> log(q)[1:10, 1:5]
           V1         V2         V3         V4         V5
1  -0.6265915 -4.6714446 -4.4833791 -1.5409610 -5.8048438
2  -3.6248309 -0.9877150 -0.8033024 -4.0219634 -0.8982287
3  -1.6434192 -5.6254270 -5.7487974 -0.7273009 -7.2508837
4  -2.1903142       -Inf       -Inf -1.4713235 -4.4292569
5  -2.4113228 -0.5256624 -0.6737870 -5.7278079 -0.5913405
6  -4.8870614       -Inf       -Inf -3.3106958 -6.6340431
7  -5.3261117       -Inf -6.1544972 -4.2220715 -7.8107129
8  -4.6803038 -4.4472730 -4.2427592       -Inf       -Inf
9  -4.9757744 -5.4457675 -4.4737512       -Inf -4.3670203
10 -4.1698733 -4.9034546       -Inf -7.3217241 -4.7028016

您可以看到，大量的-Inf值使我的计算混乱。我曾考虑过用很小的数字代替零，但是总和不再是1。我如何编写代码以构造q的替代矩阵，即1）没有零值，因此绕过了log（0）问题，并且2）仍然具有加起来为一的列，并且不更改...的基本分布行中的数据？

非常感谢！

编辑：为了提供更广泛的上下文：我需要进行对数转换，因为我正在将输出提供给计算对数似然函数。在我的应用程序中，我正在重新参数化Dirichlet分布的对数似然性，因此我没有从包中调用通用对数似然函数。

这是我的整体功能：

llikelihood = function(alpha0, beta, q, d, n) {
  llike = n*(lgamma(alpha0) - sum_a(alpha0, beta, d) + sum_b (alpha0, beta, q, d, n))
  return(llike)
}

sum_a = function(alpha0, beta, d) {
  sum_a = 0
  for (i in 1:d) {
    sum_a = sum_a + lgamma(alpha0*beta[i]) 
  }
  return(sum_a)
}

# returns the output to summation from 1 to k of (alpha0*beta[i] - 1)*log(x_i)
sum_b = function(alpha0, beta, q, d, n) { 
  # replace zero values
  sum_b = 0
  # find the log q
  logq = log(q)
  qlog = apply(logq, 1, sum)
  #  for each column, sum up the draws
  for (i in 1:d) {
    sum_b = sum_b + (alpha0*beta[i] - 1)*1/n*qlog[i]
  }

  # apply(log(q), 2, sum)
  return(sum_b)
}

这里，sum_b是我如上所述计算log（q）的位置。如您所见，我的问题是我需要摆脱零，将数据归一化，然后取对数。如何编写能有效执行的代码？我猜想它就像Laplace Smoothing一样，但是我对此不太了解，并且是R语言编程的新手。非常感谢您的评论！

Answer 1

1]您可以尝试不返回零的-Inf的其他转换，例如平方根或立方根。

2）通过将所有元素除以它们的列总和来归一化1）的结果。

set.seed(123)
X <- t(rdirichlet(4, alpha=c(1,0,2,1)))
X
           [,1]      [,2]      [,3]       [,4]
[1,] 0.03562445 0.3384606 0.5700819 0.01357789
[2,] 0.00000000 0.0000000 0.0000000 0.00000000
[3,] 0.64748450 0.2927702 0.3297736 0.88378152
[4,] 0.31689105 0.3687692 0.1001445 0.10264059

colSums(X)
# [1] 1 1 1 1

步骤1）平方根。

X2 <- sqrt(X); X2
          [,1]      [,2]      [,3]      [,4]
[1,] 0.1887444 0.5817737 0.7550377 0.1165242
[2,] 0.0000000 0.0000000 0.0000000 0.0000000
[3,] 0.8046642 0.5410824 0.5742592 0.9400966
[4,] 0.5629308 0.6072637 0.3164561 0.3203757

步骤2）标准化

X3 <- sweep(X2, 2, colSums(X2), FUN="/"); X3

          [,1]      [,2]      [,3]       [,4]
[1,] 0.1212746 0.3362621 0.4587794 0.08462201
[2,] 0.0000000 0.0000000 0.0000000 0.00000000
[3,] 0.5170236 0.3127428 0.3489340 0.68271531
[4,] 0.3617018 0.3509952 0.1922865 0.23266269

> colSums(X3)
[1] 1 1 1 1