R中线性判别分析中的分类函数答案

【问题标题】：Classification functions in linear discriminant analysis in RR中线性判别分析中的分类函数
【发布时间】：2011-08-03 12:38:33
【问题描述】：

在 R 中使用 lda() 完成线性判别分析后，是否有一种方便的方法可以为每个组提取 classification functions？

从链接，

不要将这些与判别函数混淆。分类函数可用于确定每个案例最可能属于哪个组。有多少组就有多少分类功能。每个函数都允许我们通过应用公式计算每个组的每个案例的分类分数：

Si = ci + wi1*x1 + wi2*x2 + ... + wim*xm

在这个公式中，下标i表示相应的组；下标 1, 2, ..., m 表示 m 个变量； ci 是第 i 组的常数，wij 是计算第 i 组分类分数时第 j 个变量的权重； xj 是第 j 个变量的相应情况的观察值。 Si 是最终的分类分数。

我们可以使用分类函数直接计算一些新观察的分类分数。

我可以使用教科书公式从头开始构建它们，但这需要从 lda 分析中重建许多中间步骤。有没有办法事后从 lda 对象中获取它们？

添加：

除非我仍然对 Brandon 的回答有误解（抱歉造成混淆！），否则答案似乎是否定的。想必大部分用户都可以从predict()获取他们需要的信息，predict()提供了基于lda()的分类。

【问题讨论】：

可以从MASS:::predict.lda 的代码中提取分类函数。实际上有三个不同的版本。默认方法是使用“插件”。还有一个附加术语log(priors)，我在上面没有看到。还有一个求幂步骤，但距离度量应在凸函数变换下保持其所需的属性。我认为可能需要这样一个步骤来维护作为结果一部分的 $posterior 矩阵的 rowSums ==1。

标签： r classification

【解决方案1】：

假设 x 是你的 LDA 对象：

x$terms

您可以通过查看对象的结构来获得一个峰值：

str(x)

更新：

Iris <- data.frame(rbind(iris3[,,1], iris3[,,2], iris3[,,3]),Sp = rep(c("s","c","v"), rep(50,3)))
train <- sample(1:150, 75)
table(Iris$Sp[train])
z <- lda(Sp ~ ., Iris, prior = c(1,1,1)/3, subset = train)
predict(z, Iris[-train, ])$class
str(z)
List of 10
 $ prior  : Named num [1:3] 0.333 0.333 0.333
  ..- attr(*, "names")= chr [1:3] "c" "s" "v"
 $ counts : Named int [1:3] 30 25 20
  ..- attr(*, "names")= chr [1:3] "c" "s" "v"
 $ means  : num [1:3, 1:4] 6.03 5.02 6.72 2.81 3.43 ...
  ..- attr(*, "dimnames")=List of 2
  .. ..$ : chr [1:3] "c" "s" "v"
  .. ..$ : chr [1:4] "Sepal.L." "Sepal.W." "Petal.L." "Petal.W."
 $ scaling: num [1:4, 1:2] 0.545 1.655 -1.609 -3.682 -0.443 ...
  ..- attr(*, "dimnames")=List of 2
  .. ..$ : chr [1:4] "Sepal.L." "Sepal.W." "Petal.L." "Petal.W."
  .. ..$ : chr [1:2] "LD1" "LD2"
 $ lev    : chr [1:3] "c" "s" "v"
 $ svd    : num [1:2] 33.66 2.93
 $ N      : int 75
 $ call   : language lda(formula = Sp ~ ., data = Iris, prior = c(1, 1, 1)/3, subset = train)
 $ terms  :Classes 'terms', 'formula' length 3 Sp ~ Sepal.L. + Sepal.W. + Petal.L. + Petal.W.
  .. ..- attr(*, "variables")= language list(Sp, Sepal.L., Sepal.W., Petal.L., Petal.W.)
  .. ..- attr(*, "factors")= int [1:5, 1:4] 0 1 0 0 0 0 0 1 0 0 ...
  .. .. ..- attr(*, "dimnames")=List of 2
  .. .. .. ..$ : chr [1:5] "Sp" "Sepal.L." "Sepal.W." "Petal.L." ...
  .. .. .. ..$ : chr [1:4] "Sepal.L." "Sepal.W." "Petal.L." "Petal.W."
  .. ..- attr(*, "term.labels")= chr [1:4] "Sepal.L." "Sepal.W." "Petal.L." "Petal.W."
  .. ..- attr(*, "order")= int [1:4] 1 1 1 1
  .. ..- attr(*, "intercept")= int 1
  .. ..- attr(*, "response")= int 1
  .. ..- attr(*, ".Environment")=<environment: R_GlobalEnv
  .. ..- attr(*, "predvars")= language list(Sp, Sepal.L., Sepal.W., Petal.L., Petal.W.)
  .. ..- attr(*, "dataClasses")= Named chr [1:5] "factor" "numeric" "numeric" "numeric" ...
  .. .. ..- attr(*, "names")= chr [1:5] "Sp" "Sepal.L." "Sepal.W." "Petal.L." ...
 $ xlevels: Named list()
 - attr(*, "class")= chr "lda"

【讨论】：

我认为这里缺少一些东西：我的 lda 对象中没有 terms 元素！
你在使用 MASS 包中的 lda() 吗？
是的，并且没有术语 - str() 显示：“prior”、“counts”、“means”、“scaling”、“lev”、“svd”、“N”、“call "
奇怪，我在 ?lda 底部的示例中使用了 Iris 数据。也许，如果您能提供一些示例数据？
我现在正在查看来自 R 2.12.2 的 lda 示例，并且没有提及术语。您使用的是哪个版本的 R？

【解决方案2】：

我认为你的问题有缺陷......好吧，也许没有缺陷，但至少有点误导。判别函数是指组之间的距离，因此没有与单个组相关联的函数，而是描述任何两个组质心之间距离的函数。我只是answered a more recent question 并放置了一个使用 iris 数据集计算得分函数并使用它在预测变量的二维图中标记案例的示例。在 2 组分析的情况下，一个组的函数将大于零，而另一组的函数将小于零。

【讨论】：

我没有问判别函数，我问的是分类函数。它们不是同一件事。有关说明，请参阅链接。
对。分类函数将数据投影到判别函数上，并根据阈值规则对其进行离散化。如果你想看看 MASS::predict.lda 是如何做到的，那么加载 MASS 并输入getAnywhere(predict.lda)。
lda 中的每个组都有一个分类函数，如我提供的链接以及我引用的书所示。我现在在链接页面中引用了一段话来澄清一下，这样我们至少可以确定我们在争论同一件事。 predict.lda 提供应用分类函数的结果，但不返回函数本身。这就是我要求的。我仍然不明白这个问题有什么缺陷或误导。
我正在寻找同样的东西 (stackoverflow.com/questions/30495710/…)。在 SAS 中，从 LDA 中获取非常容易。(support.sas.com/documentation/cdl/en/statug/63033/HTML/default/…)

【解决方案3】：

没有内置的方法来获取我需要的信息，所以我写了一个函数来做：

ty.lda <- function(x, groups){
  x.lda <- lda(groups ~ ., as.data.frame(x))

  gr <- length(unique(groups))   ## groups might be factors or numeric
  v <- ncol(x) ## variables
  m <- x.lda$means ## group means

  w <- array(NA, dim = c(v, v, gr))

  for(i in 1:gr){
    tmp <- scale(subset(x, groups == unique(groups)[i]), scale = FALSE)
    w[,,i] <- t(tmp) %*% tmp
  }

  W <- w[,,1]
  for(i in 2:gr)
    W <- W + w[,,i]

  V <- W/(nrow(x) - gr)
  iV <- solve(V)

  class.funs <- matrix(NA, nrow = v + 1, ncol = gr)
  colnames(class.funs) <- paste("group", 1:gr, sep=".")
  rownames(class.funs) <- c("constant", paste("var", 1:v, sep = "."))

  for(i in 1:gr) {
    class.funs[1, i] <- -0.5 * t(m[i,]) %*% iV %*% (m[i,])
    class.funs[2:(v+1) ,i] <- iV %*% (m[i,])
  }

  x.lda$class.funs <- class.funs

  return(x.lda)
}

此代码遵循勒让德和勒让德的数值生态学 (1998) 第 625 页中的公式，并与从第 626 页开始的工作示例的结果相匹配。

【讨论】：