Skip to main content
Post Undeleted by Michael E2
added 540 characters in body
Source Link
Michael E2
  • 244.7k
  • 18
  • 351
  • 774

Compiled, parallelized ordering plus a more functional approach. The switch between DistanceMatrix and Outer (for dm) may be system dependent.

kMeans2[data_, k_]ordC := 
 Module[Compile[{{x =, N[data]_Real, c1}}, 
 dmat  First@Ordering[x, cn1], 
 n = Length[data],RuntimeAttributes minpos,-> clusters{Listable},
  Parallelization -> True];

kMeans2 sum// ClearAll;
kMeans2[data_, tk_] := 0Module[{x = Developer`ToPackedArray@N[data], dm},(*Randomly 
 choose kIf[ArrayDepth[x] centriods> from1 the|| datak*Length[x] points*)< 10^5,
  c dm =(*Sqrt@*)Total[Outer[Subtract, RandomSample[x##, k];
1]^2, {3}] Do[t++;&,
   (*Assigndm each= pointDistanceMatrix;
 to the cluster];
 of theNestWhile[
 closest centriod*) With[
   dmat = DistanceMatrix[x,{clusters c];= 
       With[{minpos = ordC[dmat];ordC[dm[x, Last@#]]}, 
   clusters =    Pick[x, minpos, #] & /@ Range@k;Range@k]},
   (*Calculate the new{clusters, centriods*)
Mean /@ clusters}
 cn = Mean /@ clusters;] &,
   (*If* theinitial centriods{clusters, aremeans} the*)
 "almost" same then{{}, terminate*)RandomSample[x, k]},
   If[Norm[cNorm[Last[#1] - cnLast[#2], Infinity] <=- 0.001,
  > 0 Break[]&,
    c2,
 = cn];, {10100 (*Maximum number of iterations*)]
  ];

Example:

foo = RandomReal[10, {10^3, 3}];
(SeedRandom[0];
  km1 = kMeans[foo, 8];) // RepeatedTiming
(SeedRandom[0];
  km2 = kMeans2[foo, 8];) // RepeatedTiming
km1 == km2
(*
{clusters0.212839, cNull}]
{0.0588948, Null}
True
*)
kMeans2[data_, k_] := 
 Module[{x = N[data], c, dmat, cn, n = Length[data], minpos, clusters,
     sum, t = 0},(*Randomly choose k centriods from the data points*)
  c = RandomSample[x, k];
  Do[t++;
   (*Assign each point to the cluster of the closest centriod*)
   dmat = DistanceMatrix[x, c];
   minpos = ordC[dmat];
   clusters = Pick[x, minpos, #] & /@ Range@k;
   (*Calculate the new centriods*)
   cn = Mean /@ clusters;
   (*If the centriods are the "almost" same then terminate*)
   If[Norm[c - cn, Infinity] <= 0.001,
    Break[],
    c = cn];, {10 (*Maximum number of iterations*)}];
  {clusters, c}]

Compiled, parallelized ordering plus a more functional approach. The switch between DistanceMatrix and Outer (for dm) may be system dependent.

ordC = Compile[{{x, _Real, 1}}, 
   First@Ordering[x, 1], 
   RuntimeAttributes -> {Listable}, Parallelization -> True];

kMeans2 // ClearAll;
kMeans2[data_, k_] := Module[{x = Developer`ToPackedArray@N[data], dm}, 
  If[ArrayDepth[x] > 1 || k*Length[x] < 10^5,
   dm =(*Sqrt@*)Total[Outer[Subtract, ##, 1]^2, {3}] &,
   dm = DistanceMatrix;
   ];
  NestWhile[
   With[
     {clusters = 
       With[{minpos = ordC[dm[x, Last@#]]}, 
        Pick[x, minpos, #] & /@ Range@k]},
     {clusters, Mean /@ clusters}
     ] &,
   (* initial {clusters, means} *)
   {{}, RandomSample[x, k]},
   Norm[Last[#1] - Last[#2], Infinity] - 0.001 > 0 &,
   2,
   100 (*Maximum number of iterations*)]
  ];

Example:

foo = RandomReal[10, {10^3, 3}];
(SeedRandom[0];
  km1 = kMeans[foo, 8];) // RepeatedTiming
(SeedRandom[0];
  km2 = kMeans2[foo, 8];) // RepeatedTiming
km1 == km2
(*
{0.212839, Null}
{0.0588948, Null}
True
*)
Post Deleted by Michael E2
Source Link
Michael E2
  • 244.7k
  • 18
  • 351
  • 774

kMeans2[data_, k_] := 
 Module[{x = N[data], c, dmat, cn, n = Length[data], minpos, clusters,
    sum, t = 0},(*Randomly choose k centriods from the data points*)
  c = RandomSample[x, k];
  Do[t++;
   (*Assign each point to the cluster of the closest centriod*)
   dmat = DistanceMatrix[x, c];
   minpos = ordC[dmat];
   clusters = Pick[x, minpos, #] & /@ Range@k;
   (*Calculate the new centriods*)
   cn = Mean /@ clusters;
   (*If the centriods are the "almost" same then terminate*)
   If[Norm[c - cn, Infinity] <= 0.001,
    Break[],
    c = cn];, {10 (*Maximum number of iterations*)}];
  {clusters, c}]