在R中的线性判别分析中的分类功能

Question

在R中的线性判别分析中的分类功能

9

完成在R中使用lda()进行线性判别分析后，是否有方便的方法提取每个组的分类函数？这些不应与判别函数混淆。分类函数可用于确定每个案例最可能属于哪个组。每个组都有一个分类函数。通过应用以下公式，每个函数允许我们为每个组的每个案例计算分类分数：

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()中获取所需的信息，该函数基于lda()进行分类。

- Tyler

可以将MASS:::predict.lda中的分类函数提取出来。实际上有三个不同的版本。默认方法是使用“插件”。我没有看到上面表示的log(priors)的额外项。还有一个指数步骤，但距离度量应该在通过凸函数的转换下保持其所需的属性。我认为这样的步骤可能需要维护行总和==1，因为它是结果的一部分的$posterior matrix。 - IRTFM

3个回答

1

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

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)
}

这段代码遵循Legendre和Legendre的《数值生态学》（1998年）第625页中的公式，并与从第626页开始的示例计算结果匹配。

- Tyler

0

我认为你的问题存在一些缺陷...好吧，也许不是缺陷，但至少有些误导性。判别函数指的是组之间的距离，因此没有与单个组相关联的函数，而是描述任意两个组质心之间距离的函数。我刚刚回答了一个更近期的问题，并提供了使用鸢尾花数据集计算得分函数并将其用于标记预测变量二维图中案例的示例。在进行两组分析时，该函数将对一组大于零，对另一组小于零。

- IRTFM

我没有问判别函数，我问的是分类函数。它们不是同一回事。请参考链接了解详情。 - Tyler

好的。分类函数将数据投影到判别函数上，并基于阈值规则对其进行离散化。如果您想了解MASS::predict.lda是如何实现的，请加载MASS并键入getAnywhere(predict.lda)。 - IRTFM

在lda中，每个组都有一个分类函数，正如我提供的链接和引用的书籍所示。我现在已经包含了链接页面上的一句话来澄清，这样我们至少可以确定我们正在讨论同一件事情。predict.lda提供了应用分类函数的结果，但它不会返回函数本身。这就是我所要求的。我仍然不明白这个问题有什么缺陷或误导之处。 - Tyler

我正在寻找与此相同的东西（http://stackoverflow.com/questions/30495710/linear-discriminant-analysis）。在SAS中，从LDA获取非常容易。（http://support.sas.com/documentation/cdl/en/statug/63033/HTML/default/viewer.htm#statug_discrim_sect029.htm） - Juanchi

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- Brandon Bertelsen · Accepted Answer

假设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"