当研究星系团的sd时，有时返回NA，有时没有

Question

我有一系列的观察，我用手段聚在一起。然后，对每个簇内的标准差(sd)进行了研究，得到了最大偏差。

如果我运行几次相同的代码，有时会出现NA。

减少k_n (集群数量)会产生同样的效果，增加会更糟；na.rm=T的存在或不存在的情况不会改变。

有人能解释一下我做错了什么吗？

守则：

k_n <-11
clusters=as.data.frame(kmeans(ex,k_n, nstart=50,iter.max = 15 )$cluster)
clusters<-cbind(clusters,ex)
temp<-sapply(1:k_n, function(k){temp=subset(clusters, clusters[,1]==k)
                                sd<-sd(temp[,2], na.rm = T)
                                return(sd)})
max(temp)
temp

以下是三分的结果。正如你所看到的，第三次审判返回一个NA，其他两个没有。

​

​

以下是数据"ex"：

2000 2000 20001200 1200 1200 1600 1600 1600 1200 1200 1600 1600 1600 1200 1600 1600 1600 1200 1200 2000 1600 1600 1600 1200 1200 1200 2000 2000 1600 1200 1200 1200 2000 2000 2000 1900 1900 1900 1350 2000 2000 2000 1600 1600 16002000 1900 1901 2000 2000 2000 1200 1200 1900 1900 1900 1500 1500 1500 1900 1900 1900 1350 2000 1600 1600 1900 1900 1900 1500 1600 1900 1900 1900 1500 1500 1900 1900 1900 1500 1200 1200 2000 20002000 2000 2000 1900 1900 1900 1600 1600 1600 1200 1200 1200 2000 1200 1200 1200 2000 2000 2000 1900 1900 1900 1901 1901 19011200 1600 1600 1600 1900 1900 1900 1300 1900 1900 1900 2000 2000 2000 1200 1900 1900 1350 2000 2000 2000 1200 1200 1200 1600 1600 16002000 2000 1200 1200 1200 1600 1600 1600 1900 1900 1900 1400 1400 1400 140 2000 2000 20002000 2000 2000

Answer 1

您可以获得sd的NA，因为集群中只有一个成员。

使用您的例子：

set.seed(1)
clus = kmeans(ex,k_n, nstart=50,iter.max = 15 )

使用您的代码：

clusters=as.data.frame(clus$cluster)
clusters<-cbind(clusters,ex)
temp<-sapply(1:k_n, function(k){temp=subset(clusters, clusters[,1]==k)
                                sd<-sd(temp[,2], na.rm = T)
                                return(sd)})

> temp
 [1]  0.40191848  4.10391341  0.00000000  0.06097108  1.57499631 12.59800291
 [7]  0.00000000  0.00000000          NA  0.00000000  0.00000000

其中一个是NA，如果您查看哪个集群给您NA，这个集群只有一个成员：

clusters[clusters[,1]==which(is.na(temp)),]
   clus$cluster   ex
69            9 1960

如果我们看看你的数据：

table(ex)
ex
   0  801 1100 1200 1350 1400 1401 1500 1560 1600 1800 1900 1901 1920 1960 2000 
   2    2    2  123   36   21    5   22    1   94    4  147   18    1    1  268 
2001 
   1

我认为，如果你增加k，你可能最终会有一个簇，只有一个成员，允许收敛。

我可以建议的一种方法是增加开始的次数：

STARTS = seq(50,500,by=50)
# we test over 50 reps, how many single clusters we get
n_equal_one = sapply(STARTS,function(S){
replicate(50,sum(kmeans(ex,k_n, nstart=S,iter.max = 15 )$size==1))
})

plot(STARTS,colMeans(n_equal_one),ylab="Average proportion of singleton cluster")

​

​

因此，如果您尝试nstart = 400或500，您将避免单例(带有n=1的集群)，但是如果数据变得更加稀疏，这可能是不可避免的。

dput(ex)

c(1400, 1400, 2000, 2000, 2000, 2001, 1400, 2000, 2000, 2000, 
2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 1400, 1400, 1401, 
2000, 2000, 2000, 2000, 1401, 1401, 2000, 2000, 2000, 2000, 2000, 
2000, 2000, 2000, 801, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 
2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 
2000, 1400, 1400, 2000, 2000, 2000, 1400, 1400, 2000, 2000, 2000, 
1200, 1600, 1500, 1960, 1350, 1900, 1900, 1900, 1350, 2000, 2000, 
2000, 2000, 1200, 1200, 1200, 1600, 1600, 1600, 1600, 1600, 1600, 
1600, 1600, 1600, 1200, 1200, 1200, 1200, 1200, 1200, 1200, 1200, 
1200, 2000, 2000, 2000, 2000, 2000, 2000, 1900, 1350, 1900, 1900, 
1350, 1350, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 1200, 1200, 
1200, 1200, 1200, 1200, 1200, 1600, 1600, 1600, 1600, 1600, 1600, 
1600, 1600, 1200, 1200, 1200, 1200, 1200, 1200, 1200, 2000, 2000, 
2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 1920, 1350, 1350, 
1900, 1900, 1900, 1900, 1900, 1350, 2000, 2000, 2000, 2000, 2000, 
2000, 2000, 1200, 1200, 1200, 1200, 1200, 1600, 1600, 1600, 1600, 
1600, 1600, 1600, 1600, 1600, 1600, 1200, 1200, 1200, 1200, 1200, 
1200, 1200, 1200, 1200, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 
1900, 1350, 1900, 1900, 1900, 1900, 1350, 2000, 2000, 2000, 2000, 
2000, 2000, 2000, 2000, 1200, 1200, 1600, 1600, 1600, 1600, 1600, 
1600, 1600, 1600, 1600, 1200, 1200, 1200, 1200, 1200, 1200, 1200, 
1200, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 
1900, 1901, 2000, 2000, 2000, 2000, 1200, 1200, 1200, 1200, 1200, 
1200, 1200, 1200, 1200, 1600, 1600, 1600, 1600, 1600, 1600, 1600, 
1600, 1600, 1600, 1200, 1200, 1200, 1200, 1200, 1200, 1200, 1200, 
1200, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 1900, 
1900, 1900, 1900, 1900, 1900, 1900, 1900, 1900, 1900, 1900, 1900, 
1900, 1500, 1500, 1500, 1500, 2000, 1350, 1900, 1900, 1900, 1350, 
2000, 2000, 1900, 1900, 1900, 1900, 1900, 1900, 1900, 1900, 1900, 
1900, 1900, 1900, 1900, 1900, 1900, 1900, 1900, 1500, 1500, 1400, 
1350, 1350, 1900, 1900, 1900, 1350, 2000, 2000, 2000, 2000, 2000, 
2000, 1200, 1200, 1200, 1200, 1200, 1600, 1600, 1600, 1600, 1600, 
1200, 1200, 1200, 1200, 1200, 1200, 1200, 1200, 2000, 2000, 2000, 
2000, 2000, 2000, 1900, 1900, 1900, 1900, 1900, 1900, 1900, 1900, 
1900, 1900, 1900, 1900, 1900, 1500, 1500, 1500, 1500, 0, 0, 1900, 
1900, 1900, 1900, 1900, 1900, 1900, 1900, 1900, 1900, 1900, 1900, 
1900, 1900, 1500, 1500, 1500, 2000, 2000, 2000, 2000, 2000, 1200, 
1200, 1200, 1200, 1200, 1200, 1600, 1600, 1600, 1600, 1600, 1600, 
1600, 1600, 1600, 1200, 1200, 1200, 1200, 1200, 1200, 1200, 1200, 
2000, 2000, 2000, 2000, 2000, 1900, 1900, 1900, 1900, 1900, 1900, 
1500, 1500, 1500, 1400, 1901, 1901, 1901, 1901, 2000, 1350, 1900, 
1900, 1900, 1350, 1350, 2000, 2000, 2000, 2000, 1500, 1560, 1900, 
1900, 1900, 1900, 1900, 1900, 1400, 2000, 1350, 1900, 1900, 1900, 
1350, 2000, 2000, 2000, 1200, 1200, 1600, 1600, 1600, 1600, 1600, 
1600, 1600, 1600, 1200, 1200, 1200, 1200, 1200, 1200, 1200, 2000, 
2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 1400, 1900, 
1900, 1900, 1900, 1900, 1900, 1900, 1900, 1900, 1900, 1900, 1900, 
1900, 1900, 1900, 1900, 1900, 1900, 1500, 1500, 1500, 1500, 2000, 
1900, 1900, 1900, 1900, 1900, 1900, 1900, 1900, 2000, 1350, 1900, 
1900, 1900, 1350, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 
2000, 2000, 1200, 1200, 1200, 1200, 1200, 1200, 1200, 1200, 1600, 
1600, 1600, 1600, 1600, 1600, 1600, 1600, 1600, 1600, 1600, 1600, 
1600, 1900, 1900, 1350, 1350, 1350, 1350, 1350, 1900, 2000, 1350, 
1900, 1900, 1900, 1900, 1900, 1900, 1350, 1350, 2000, 2000, 2000, 
2000, 2000, 2000, 2000, 2000, 2000, 2000, 1200, 1200, 1200, 1200, 
1200, 1200, 1200, 1200, 1200, 1600, 1600, 1600, 1600, 1600, 1600, 
1600, 1600, 1600, 1600, 1600, 1900, 1900, 801, 1400, 1400, 1400, 
1400, 1901, 1901, 1900, 1901, 1901, 1900, 1901, 1901, 1901, 1900, 
1900, 1900, 1900, 1901, 1900, 1900, 1901, 1901, 1900, 1901, 1901, 
1901, 1350, 1350, 1350, 1350, 1350, 1350, 2000, 2000, 2000, 2000, 
2000, 2000, 2000, 2000, 2000, 1400, 1400, 1401, 1401, 2000, 2000, 
2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 1100, 1100, 
2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 
2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 
2000, 2000, 2000, 2000, 1800, 2000, 1800, 1800, 1400, 1200, 1400, 
1600, 1800, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 
2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 
2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000
)