Unüberwachte Klassifizierung mit km in R.

Ich habe eine Zeitreihe von Satellitenbildern (5 Bänder) und möchte sie nach km in R klassifizieren. Mein Skript funktioniert einwandfrei (durchlaufen Sie meine Bilder, konvertieren Sie die Bilder in data.frame, gruppieren Sie sie und konvertieren Sie sie zurück in a Raster):

for (n in files) {
image <- stack(n)    
image <- clip(image,subset)

###classify raster
image.df <- as.data.frame(image)  
cluster.image <- kmeans(na.omit(image.df), 10, iter.max = 10, nstart = 25) ### kmeans, with 10 clusters

#add back NAs using the NAs in band 1 (identic NA positions in all bands), see http://stackoverflow.com/questions/12006366/add-back-nas-after-removing-them/12006502#12006502
image.df.factor <- rep(NA, length(image.df[,1]))
image.df.factor[!is.na(image.df[,1])] <- cluster.image$cluster

#create raster output
clusters <- raster(image)   ## create an empty raster with same extent than "image"  
clusters <- setValues(clusters, image.df.factor) ## fill the empty raster with the class results  
plot(clusters)
}

Mein Problem ist: Ich kann die Klassifizierungsergebnisse nicht miteinander vergleichen, da sich die Clusterzuweisungen von Bild zu Bild unterscheiden. Zum Beispiel befindet sich "Wasser" im ersten Bildcluster Nummer 1, in den nächsten 2 und in der dritten 10, was es unmöglich macht, die Wasserergebnisse zwischen den Daten zu vergleichen.

Wie kann ich die Clusterzuordnung korrigieren?

Kann ich für alle Bilder einen festen Startpunkt angeben (in der Hoffnung, dass Wasser immer zuerst erkannt und somit als 1 klassifiziert wird)?

Und wenn ja, wie?

Ich denke, Sie können nicht ... Sie müssen zuerst jede Klasse beschriften, um sie zu vergleichen. Kmean klassifiziert unbeaufsichtigt, ohne vorherige Informationen, und kann daher keine Klassen definieren.

Wenn Sie eine Referenzschicht haben, können Sie eine Kennzeichnung mit Stimmenmehrheit vornehmen. Hier ist ein wesentlich effizienterer Code für die Mehrheitsentscheidung als die Verwendung der Paketfunktion "Raster" zonal:

require (data.table)
fun <- match.fun(modal)
vals <- getValues(ref) 
zones <- round(getValues(class_file), digits = 0) 
rDT <- data.table(vals, z=zones) 
setkey(rDT, z) 
zr<-rDT[, lapply(.SD, modal,na.rm=T), by=z]

Wo refist class_fileIhre Rasterklassenreferenzdatei , ist Ihr kmeans-Ergebnis.

zr gibt Ihnen in der ersten Spalte die 'Zonen'-Nummer und in der zweiten Spalte die Bezeichnung für die Klasse.

Um das Clustering auf einem Image-Stack zu implementieren, tun Sie dies nicht bandweise, sondern gleichzeitig auf dem gesamten Image-Stack. Andernfalls macht die Statistik, wie von @nmatton hervorgehoben, wenig Sinn.

Ich stimme jedoch nicht zu, dass dies nicht möglich ist, sondern nur speicherintensiv. Bei realen Satellitendaten ist dies ein großes Problem und bei hochauflösenden Daten möglicherweise unmöglich. Sie können es jedoch im Speicher verarbeiten, indem Sie Ihre Raster in ein einzelnes Objekt zwingen, das an eine Clustering-Funktion übergeben werden kann. Sie müssen NA-Werte über Raster hinweg verfolgen, da diese während des Clusters entfernt werden und Sie die Positionen im Raster kennen müssen, damit Sie die Clusterwerte den richtigen Zellen zuweisen können.

Wir können hier einen Ansatz durchgehen. Fügen wir die erforderlichen Bibliotheken und einige Beispieldaten hinzu (das RGB R-Logo, mit dem wir 3 Bänder zum Arbeiten erhalten).

library(raster)
library(cluster)
r <- stack(system.file("external/rlogo.grd", package="raster")) 
  plot(r)

Erstens können wir unser Multiband-Raster-Stack-Objekt mithilfe von getValues ​​in einen data.frame zwingen. Beachten Sie, dass ich in Zeile 1, Spalte 3 einen NA-Wert hinzufüge, damit ich veranschaulichen kann, wie mit keinen Daten umgegangen wird.

r.vals <- getValues(r[[1:3]])
  r.vals[1,][3] <- NA

Hier können wir zur Sache kommen und einen Zellenindex der Nicht-NA-Werte erstellen, der zum Zuweisen der Clusterergebnisse verwendet wird.

idx <- 1:ncell(r)
idx <- idx[-unique(which(is.na(r.vals), arr.ind=TRUE)[,1])]  

Nun erstellen wir ein Clusterobjekt aus den 3-Band-RGB-Werten mit k = 4. Ich verwende die Clara K-Medoids-Methode, weil sie mit großen Datenmengen gut und mit ungeraden Verteilungen besser ist. Es ist K-Means sehr ähnlich.

clus <- cluster::clara(na.omit(scale(r.vals)), k=4)

Der Einfachheit halber können wir ein leeres Raster erstellen, indem wir eines der Rasterbänder aus unserem ursprünglichen Rasterstapelobjekt ziehen und ihm NA-Werte zuweisen.

r.clust <- r[[1]]
r.clust[] <- NA

Schließlich weisen wir mithilfe des Index die Clusterwerte der entsprechenden Zelle im leeren Raster zu und zeichnen die Ergebnisse auf.

r.clust[idx] <- clus$clustering
plot(r.clust) 

Für große Raster sollten Sie sich das BigMemory-Paket ansehen, das Matrizen auf die Festplatte schreibt und Blöcke bearbeitet, und es steht eine k-means-Funktion zur Verfügung. Beachten Sie auch, dass dies nicht genau das ist, wofür R entwickelt wurde, und dass eine Bildverarbeitung oder GIS-Software möglicherweise besser geeignet ist. Ich weiß, dass SAGA und die Orfeo-Toolbox beide freie Software sind, für die K-Means-Clustering für Image-Stacks verfügbar ist. Es gibt sogar eine RSAGA-Bibliothek, mit der die Software von R aus aufgerufen werden kann.