Support Vector Machine Estimation for Japanese Flag(A circle inside a rectangle)

Support Vector Machine (SVM) can sometimes surprise with its capability. I was thinking to test it a nonlinear example and suddenly taught,why do not i model Japanese flag. I can calculate points of rectangle and then points of circle inside. Then I taught why do i calculate, i better to draw picture and extract rgb values from it. And I can even extract my data from picture. For japanese flag I used 3 colors. White and Red : as in Japanese Flag Green : part that i cut from flag. I want SVM to predict this.

I extract values and dump into 2 separate csv files.

 
public static ArrayList<ArrayList<Object>> getImageData() {

		BufferedImage img;
		ArrayList<ArrayList<Object>> result = new ArrayList<>();
		try {
			img = ImageIO.read(new File(IMG));

			int widthSample = 1;
			int heightSample = 2;

			for (int j = 0; j < img.getHeight(); j++) {				
				if (j % heightSample != 0) {
					continue;
				}
				for (int i = 0; i < img.getWidth(); i++) {

					if (i % widthSample != 0) {
						continue;
					}

					String color = getPixelData(img, i, j);

					ArrayList<Object> res = new ArrayList<>();
					res.add(i);
					res.add(j);
					res.add(color);

					result.add(res);

				}
			}

		} catch (IOException e) {
			e.printStackTrace();
		}

		return result;
	}

	private static String getPixelData(BufferedImage img, int x, int y) {
		int argb = img.getRGB(x, y);

		int rgb[] = new int[] { (argb >> 16) & 0xff, // red
				(argb >> 8) & 0xff, // green
				(argb) & 0xff // blue
		};

		String color = colorUtils.getColorNameFromRgb(rgb[0], rgb[1], rgb[2]);

		return color;
	}

With this logic, you can generate some shapes in a plane and make some part of image as green. Green part will be your test data.(Place you want svm to predict) For example you can paint 1/4 of a circle or rectangle to green to see the result. In next post i will try stars and top of bottle.

maindata <- read.csv(file="E:/REDA/japan_flag_data1.csv", header=FALSE, sep=",")


test  <- read.csv(file="E:/REDA/japan_flag_data1_test.csv", header=FALSE, sep=",")

traindata = data.frame(
  x=maindata$V1, y=maindata$V2 ,colpart = as.factor( maindata$V3 )
)

testdf = data.frame(
  x=test$V1, y=test$V2 ,colpart = as.factor( test$V3 )
)


plot(traindata$x, traindata$y , col =  traindata$colpart   )

fitsvm =  svm(colpart ~ ., data=traindata)

plot(fitsvm, traindata, y~x , col=c("red","white")) 




predictdata = predict(fitsvm, newdata = testdf)



points(x=testdf$x, y = testdf$y , pch = 19,col = c("blue", "green")[as.numeric(predictdata)] ,cex=2.6)

The blue and green points are predictions from missing part. Blue ones which are nearer to center of circle are predicted correct. Green one are not predicted correct.We can tune parameters to estimate missing values.

Residual vs Fitted for investment

While I was checking graphs i can draw from a regression model I realized
one graph is extremely useful for investment.
Think you have(actually we have) thousand of used car prices. And you want to buy a car
with optimum price for investment.(In hope you will sell later)

Best Car to buy Price = Max(Expected Price by regression - Real price )

Think our regression line expects a car price to be 30K but actual advertisement price is 20K.
There are 2 possibilities.
1)Car is damaged.
2)Car owner needs urgent money and selling his car with a very low price.

If car price is 40K. I can not find a logical explanation for this. Some people
are trying to sell their used cars with more price than an-unused one. Probably
they spend for some amenities which they think so valuable.

import org.apache.spark.ml.feature.VectorAssembler
import org.apache.spark.ml.linalg.Vectors

val dataset = spark.createDataFrame(
  Seq(
  (20000,2011,30000.0),
  (120000,2014,20000.0),
  (60000,2015,25000.0) ,
  (20000,2011,32000.0),
  (120000,2014,21000.0),
  (60000,2015,45000.0)   
  
  )
).toDF("km", "year", "label")

val assembler = new VectorAssembler()
  .setInputCols(Array("km", "year"))
  .setOutputCol("features")

val output = assembler.transform(dataset)
output.select("features", "label").show(false)

import org.apache.spark.ml.regression.LinearRegression

val lr = new LinearRegression()
  .setMaxIter(10)
  .setRegParam(0.3)
  .setElasticNetParam(0.8)


val lrModel = lr.fit(output)

display(lrModel, output, "fittedVsResiduals")

println(s"Coefficients: ${lrModel.coefficients} Intercept: ${lrModel.intercept}")

​



val trainingSummary = lrModel.summary

println(s"numIterations: ${trainingSummary.totalIterations}")

println(s"objectiveHistory: [${trainingSummary.objectiveHistory.mkString(",")}]")

trainingSummary.residuals.show()

println(s"RMSE: ${trainingSummary.rootMeanSquaredError}")

println(s"r2: ${trainingSummary.r2}")

Coefficients: [-0.19283786035452538,2928.112739446878] Intercept: -5853577.79139608
numIterations: 8
objectiveHistory: [0.5,0.446369757257132,0.352850077757605,0.272318835721877,0.26365142412164966,0.23726105027025182,0.23725993458647637,0.23725993457463043]
+-------------------+
|          residuals|
+-------------------+
|-1000.1704245014116|
|-500.72260739002377|
| -9999.106968107633|
|  999.8295754985884|
| 499.27739260997623|
| 10000.893031892367|
+-------------------+

Residual = Observed value - Predicted value

We must find the ones with most negative residual.(Much cheaper than expected).
Databricks graph has bad resolution for few points so i wrote R version also.

library(lattice) 
mydata2 = data.frame(
  year = c(2011.0,2012.0,2014.0,2015.0),
  km10000 = c(6.0,7.0,10.0,3.0),
  price1000 = c(200.0,250.0,300.0,400.0)
)



res2.mod1 = lm(price1000 ~  km10000 + year , data = mydata2)
summary(res.mod1)
fitted(res2.mod1)
xyplot(resid(res2.mod1) ~ fitted(res2.mod1),
       xlab = "Fitted Values",
       ylab = "Residuals",
       main = "Car price based on year and km ",
       par.settings = simpleTheme(col=c("blue","red"),
                                  pch=c(10,3,11), cex=3, lwd=2),
       
       panel = function(x, y, ...)
       {
         panel.grid(h = -1, v = -1)
         panel.abline(h = 0)
         panel.xyplot(x, y, ...)
       }
)       

> fitted(res2.mod1)
       1        2        3        4 
206.1722 240.9232 302.5415 400.3631 
> resid(res2.mod1)
         1          2          3          4 
-6.1721992  9.0767635 -2.5414938 -0.3630705 
> 

Impute Outliers In Spark

Imputing data in R is so easy.

fun <- function(x){
    quantiles <- quantile( x, c(.05, .95 ) )
    x[ x < quantiles[1] ] <- quantiles[1]
    x[ x > quantiles[2] ] <- quantiles[2]
    x
}
fun( yourdata )

I tried to write my functions for Spark. Below you will find methods for replacing values
Q1 is 1st quantile and
Q3 is 3rd quantile
left outliers < Q1 – 1.5×IQR
right outliers > Q3 + 1.5×IQR
For a dataset I obtain below values for these results.

org.apache.spark.util.StatCounter = (count: 13, mean: 131.538462, stdev: 65.441242, max: 320.000000, min: 1.000000) 
rddMin: Double = 1.0 
rddMax: Double = 320.0 
rddMean: Double = 131.53846153846155 
quantiles: Array[Double] = Array(117.0, 150.0) 
Q1: Double = 117.0 
Q3: Double = 150.0 
IQR: Double = 33.0 
lowerRange: Double = 67.5 
upperRange: Double = 199.5

So according to our data set and our logic we can make below changes.(Of course there are more)
left outliers -> Q1
left outliers -> mean
right outliers -> Q3
right outliers -> mean
You can also want to impute values less than Q1 to Q1.

import org.apache.spark.sql.functions._

val imputeLessToTarget: (Double,Double,Double ) => Double = (arg: Double,treshold: Double,target: Double) => { if (arg < treshold) target else arg}
val imputeLessToTargetUDF = udf(imputeLessToTarget)

val imputeMoreToTarget: (Double,Double,Double ) => Double = (arg: Double,treshold: Double,target: Double) => { if (arg > treshold) target else arg}
val imputeMoreToTargetUDF = udf(imputeMoreToTarget)


//Change values less than lowerRange of outliers to Q1 **** Capping
val df5 = df.withColumn("replaced", imputeLessToTargetUDF(col("original") ,lit(lowerRange),lit(Q1) ))
println( "df5 Change values less than lowerRange of outliers to Q1 **** Capping")
df5.show()
//Change values less than lowerRange of outliers to mean
val df6 = df.withColumn("replaced", imputeLessToTargetUDF(col("original") ,lit(lowerRange),lit(rddMean) ))
println( "df6 Change values less than lowerRange of outliers to mean")
df6.show()
//Change values less than Q1 to mean
val df7 = df.withColumn("replaced", imputeLessToTargetUDF(col("original") ,lit(Q1),lit(rddMean) ))
println( "df7 Change values less than Q1 to mean")
df7.show()
//Change values more than upperRange of outliers to Q3 **** Capping
val df8 = df.withColumn("replaced", imputeMoreToTargetUDF(col("original") ,lit(upperRange),lit(Q3) ))
println( "df8 Change values more than upperRange of outliers to Q3 **** Capping")
df8.show()
//Change values more than upperRange of outliers to mean
val df9 = df.withColumn("replaced", imputeMoreToTargetUDF(col("original") ,lit(upperRange),lit(rddMean) ))
println( "df9 Change values more than upperRange of outliers to mean")
df9.show()

At below you can see output of generating new columns. These type of functions are my toolset for dealing with new data.

df5 Change values less than lowerRange of outliers to Q1 **** Capping 
1.0| 117.0
110.0| 110.0
111.0| 111.0
112.0| 112.0
117.0| 117.0
118.0| 118.0
120.0| 120.0
122.0| 122.0
129.0| 129.0
140.0| 140.0
150.0| 150.0
160.0| 160.0
320.0| 320.0
df6 Change values less than lowerRange of outliers to mean 
1.0|131.53846153846155
110.0| 110.0
111.0| 111.0
112.0| 112.0
117.0| 117.0
118.0| 118.0
120.0| 120.0
122.0| 122.0
129.0| 129.0
140.0| 140.0
150.0| 150.0
160.0| 160.0
320.0| 320.0
df7 Change values less than Q1 to mean 
1.0|131.53846153846155
110.0|131.53846153846155
111.0|131.53846153846155
112.0|131.53846153846155
117.0| 117.0
118.0| 118.0
120.0| 120.0
122.0| 122.0
129.0| 129.0
140.0| 140.0
150.0| 150.0
160.0| 160.0
320.0| 320.0
df8 Change values more than upperRange of outliers to Q3 **** Capping
1.0| 150.0
110.0| 150.0
111.0| 150.0
112.0| 150.0
117.0| 150.0
118.0| 150.0
120.0| 150.0
122.0| 150.0
129.0| 150.0
140.0| 150.0
150.0| 150.0
160.0| 150.0
320.0| 320.0
df9 Change values more than upperRange of outliers to mean 
1.0|131.53846153846155
110.0|131.53846153846155
111.0|131.53846153846155
112.0|131.53846153846155
117.0|131.53846153846155
118.0|131.53846153846155
120.0|131.53846153846155
122.0|131.53846153846155
129.0|131.53846153846155
140.0|131.53846153846155
150.0|131.53846153846155
160.0|131.53846153846155
320.0| 320.0

R Impute Dataframe( Replace outliers )

There are some methods over internet for imputing outliers.
I here give a sample I use,which i combined methods i found.
My intention is selectively applying imputing to numeric columns.
When
Q1 is 1st quantile and
Q3 is 3rd quantile
Below ranges are outliers by definition.
below Q1 – 1.5×IQR or above Q3 + 1.5×IQR


Below code replaces value as x < min or x > max with mean value. According to your
data median could be a better choice.

numcol <- c(1,3,40,50,600)
numcol2 <- c(2,420,400,500,600)
charcol <- c("a","a","b","b","a")


df <- data.frame(a=numcol,b=charcol,c=numcol2)
#select numeric columns to change
columnsToChange <- c("a","c")

df
for(i in columnsToChange){
  Q1 <- quantile(df[,i],0.75, na.rm=TRUE) 
  max <- Q1 + (IQR(df[,i], na.rm=TRUE) * 1.5 )
  
  Q3 <- quantile(df[,i],0.25, na.rm=TRUE)
  min <- Q3 - (IQR(df[,i], na.rm=TRUE) * 1.5 )
  
  message(sprintf("min ,  max  mean  %s %s mean of column %s \n", min,max ,mean(df[,i] )) )
  
  indexesstochange <- which(df[,i] < min | df[,i] > max)
  
  message(sprintf("indexes to change %s \n", indexesstochange ))
 
  df[,i][indexesstochange] <- mean(df[,i])
}
df

It produces the output below. For column a outlier is max value at index 5. For column c outlier is at min value at index 1.

1   1 a   2
2   3 a 420
3  40 b 400
4  50 b 500
5 600 a 600

min ,  max  mean  -67.5 120.5 mean of column 138.8 

indexes to change 5 

min ,  max  mean  250 650 mean of column 384.4 

indexes to change 1 


      a   b     c
1   1.0 a 384.4
2   3.0 a 420.0
3  40.0 b 400.0
4  50.0 b 500.0
5 138.8 a 600.0

Effect of Outliers with R

***
Spark or ( R ) samples are for big files which contains thousands of lines.
Also you do not know data and can not play with it.
I put here simplest data set for spark mllib so that one can play and understand what metrics
are effected from which parameters.
It is not for seniors but perfect for beginners of who need to calibrate parameters with simple sets.
***
R has a built in dataset as cars.

carssub
   speed dist
1      4    2
2      4   10
3      7    4
4      7   22
5      8   16
6      9   10
7     10   18
8     10   26
9     10   34
10    11   17
11    11   28
12    12   14
13    12   20
14    12   24
15    12   28
16    13   26
17    13   34
18    13   34
19    13   46
20    14   26
21    14   36
22    14   60
23    14   80
24    15   20
25    15   26
26    15   54
27    16   32
28    16   40
29    17   32
30    17   40

We first take a subset of data.
Then we add outliers. (1 and 5 outliers)
Then we see the effect of outliers.
lm : for fitting linear models.
abline : add line to plot

If you look at picture you will see how a line fits to data when no outlier.
When we add only 1 outlier it changes a lot. When 5 is added it gets much more worser.

carssub <- cars[1:30, ]  # original data

carssub <- cars[1:30, ]  # original data

cars_outliers1 <- data.frame(speed=c(20), dist=c( 218))  # introduce outliers.

cars_outliers5 <- data.frame(speed=c(19,19,20,20,20), dist=c(190, 186, 210, 220, 218))  # introduce outliers.

cars_outliers10 <- data.frame(speed=c(19,19,20,20,20,21,22,23,24,25), dist=c(190, 186, 210, 220, 218,220,224,230,235,240))

cars_outliers15 <- data.frame(speed=c(19,19,20,20,20,21,22,23,24,25,26,27,28,29,30), dist=c(190, 186, 210, 220, 218,220,224,230,235,240,244,245,248,250,252))

cars_outliers20 <- data.frame(speed=c(19,19,20,20,20,21,22,23,24,25,26,27,28,29,30,31,32,33,34,35), dist=c(190, 186, 210, 220, 218,220,224,230,235,240,244,245,248,250,252,254,256,258,260,262))

cars_total_1 <- rbind(carssub, cars_outliers1)  # data with outliers.
cars_total_5 <- rbind(carssub, cars_outliers5)  # data with outliers.
cars_total_10 <- rbind(carssub, cars_outliers10)  # data with outliers.
cars_total_15 <- rbind(carssub, cars_outliers15)  # data with outliers.
cars_total_20 <- rbind(carssub, cars_outliers20)  # data with outliers.


par(mfrow=c(2, 3))
plot(carssub$speed, carssub$dist, xlim=c(0, 40), ylim=c(0, 300), main="Pure data", xlab="speed", ylab="dist", pch="*", col="red", cex=2)
abline(lm(dist ~ speed, data=carssub), col="blue", lwd=3, lty=2)

#aykiri gozlemsiz model
plot(cars_total_1$speed, cars_total_1$dist, xlim=c(0, 40), ylim=c(0, 300), main="1 outlier added", xlab="speed", ylab="dist", pch="*", col="red", cex=2)
lm1 <- lm(dist ~ speed, data=cars_total_1)
abline(lm1, col="blue", lwd=3, lty=2)
summary(lm1)

plot(cars_total_5$speed, cars_total_5$dist, xlim=c(0, 40), ylim=c(0, 300), main="5 outliers added", xlab="speed", ylab="dist", pch="*", col="red", cex=2)
lm2 <- lm(dist ~ speed, data=cars_total_5)
abline(lm2, col="blue", lwd=3, lty=2)
summary(lm2)

plot(cars_total_10$speed, cars_total_10$dist, xlim=c(0, 40), ylim=c(0, 300), main="10 outlier added", xlab="speed", ylab="dist", pch="*", col="red", cex=2)
lm10 <- lm(dist ~ speed, data=cars_total_10)
abline(lm10, col="blue", lwd=3, lty=2)
summary(lm10)


plot(cars_total_15$speed, cars_total_15$dist, xlim=c(0, 40), ylim=c(0, 300), main="15 outlier added", xlab="speed", ylab="dist", pch="*", col="red", cex=2)
lm15 <- lm(dist ~ speed, data=cars_total_15)
abline(lm15, col="blue", lwd=3, lty=2)
summary(lm15)

plot(cars_total_20$speed, cars_total_20$dist, xlim=c(0, 40), ylim=c(0, 300), main="20 outlier added", xlab="speed", ylab="dist", pch="*", col="red", cex=2)
lm20 <- lm(dist ~ speed, data=cars_total_20)
abline(lm20, col="blue", lwd=3, lty=2)
summary(lm20)

no outliers: dist = speed * 2.9 - 6.8
1 outlier: dist = speed * 6.1 -40
5 outliers: dist = speed * 11.4 - 95
You can see that slope and intercept are getting worser values up to 10 but after that it is changing shape and error function is decreasing. And outliers are also becoming a large group and they are not outliers..

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  -6.8446     8.7420  -0.783 0.440223    
speed         2.9730     0.7046   4.219 0.000233 ***

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  -40.015     19.271  -2.076 0.046817 *  
speed          6.131      1.515   4.048 0.000351 ***

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  -95.203     24.507  -3.885 0.000466 ***
speed         11.437      1.793   6.379 3.18e-07 ***

Call:
lm(formula = dist ~ speed, data = cars_total_1)

Residuals:
    Min      1Q  Median      3Q     Max 
-32.210 -15.883  -5.555   6.641 135.398 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  -40.015     19.271  -2.076 0.046817 *  
speed          6.131      1.515   4.048 0.000351 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 30.63 on 29 degrees of freedom
Multiple R-squared:  0.361, Adjusted R-squared:  0.339 
F-statistic: 16.38 on 1 and 29 DF,  p-value: 0.0003514


Call:
lm(formula = dist ~ speed, data = cars_total_5)

Residuals:
    Min      1Q  Median      3Q     Max 
-67.220 -27.755  -7.473  19.428  86.471 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  -95.203     24.507  -3.885 0.000466 ***
speed         11.437      1.793   6.379 3.18e-07 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 43.88 on 33 degrees of freedom
Multiple R-squared:  0.5522, Adjusted R-squared:  0.5386 
F-statistic: 40.69 on 1 and 33 DF,  p-value: 3.175e-07



Call:
lm(formula = dist ~ speed, data = cars_total_10)

Residuals:
    Min      1Q  Median      3Q     Max 
-81.855 -30.503  -0.082  34.346  77.691 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) -123.552     20.715  -5.964 6.37e-07 ***
speed         13.965      1.366  10.221 1.85e-12 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 44.15 on 38 degrees of freedom
Multiple R-squared:  0.7333, Adjusted R-squared:  0.7263 
F-statistic: 104.5 on 1 and 38 DF,  p-value: 1.852e-12



Call:
lm(formula = dist ~ speed, data = cars_total_15)

Residuals:
    Min      1Q  Median      3Q     Max 
-78.833 -30.297  -3.225  31.292  71.650 

Coefficients:
             Estimate Std. Error t value Pr(>|t|)    
(Intercept) -114.7213    16.5388  -6.936 1.59e-08 ***
speed         13.2679     0.9684  13.701  < 2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 42.11 on 43 degrees of freedom
Multiple R-squared:  0.8136, Adjusted R-squared:  0.8093 
F-statistic: 187.7 on 1 and 43 DF,  p-value: < 2.2e-16




Call:
lm(formula = dist ~ speed, data = cars_total_20)

Residuals:
   Min     1Q Median     3Q    Max 
-73.15 -31.29  -6.27  33.74  79.83 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) -93.3047    14.5507  -6.412 5.86e-08 ***
speed        11.6738     0.7548  15.466  < 2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 42.92 on 48 degrees of freedom
Multiple R-squared:  0.8329, Adjusted R-squared:  0.8294 
F-statistic: 239.2 on 1 and 48 DF,  p-value: < 2.2e-16

Spark and R K-Means classification

***
Spark samples are for big files which contains thousands of lines.
Also you do not know data and can not play with it.
I put here simplest data set for spark mllib so that one can play and understand what metrics
are effected from which parameters.
It is not for seniors but perfect for beginners of who need to calibrate parameters with simple sets.
***
Below code is from sample Spark documentation. I changed Rdd so that one can play and understand
how data is distributed.




Here you can play with values and observe distribution of clusters.
Always print cluster centers. It will give you a clue for large datasets.

You can easily play with dataset and number of demanded clusters to get an idea of how
K-means work.

import org.apache.spark.mllib.clustering.{KMeans, KMeansModel}
import org.apache.spark.mllib.linalg.Vectors

val parsedData = sc.parallelize(Seq(
  ( Vectors.dense(1.0, 1.0)),
  ( Vectors.dense(40.0, 40.0)),
  ( Vectors.dense(60.0, 60.0)),
  ( Vectors.dense(101.0, 101.1))
))

// Cluster the data into two classes using KMeans
val numClusters = 2
val numIterations = 20
val clusters = KMeans.train(parsedData, numClusters, numIterations)

// Evaluate clustering by computing Within Set Sum of Squared Errors
val WSSSE = clusters.computeCost(parsedData)
println("Within Set Sum of Squared Errors = " + WSSSE )
val clusterCenters = clusters.clusterCenters.map(_.toArray)
println("The Cluster Centers are = " + clusterCenters)
parsedData.collect().map( s=> println( "cluster "+clusters.predict(s) +" "+s.toString() ) )

Result

Within Set Sum of Squared Errors = 3874.8066666666673
clusterCenters: Array[Array[Double]] = Array(Array(67.0, 67.03333333333333), Array(1.0, 1.0))
cluster 1 [1.0,1.0] cluster 0 [40.0,40.0] cluster 0 [60.0,60.0] cluster 0 [101.0,101.1]

Same Code In R

pointx = c(1,2, 50, 51) 
pointy = c(1,2,50,51) 
df = data.frame(pointx, pointy)
library(ggplot2)
ggplot(df, aes(pointx, pointy)) + geom_point()
myCluster <- kmeans(df, 3, nstart = 20)
myCluster$centers
myCluster$clus <- as.factor(myCluster$cluster)
ggplot(df, aes(pointx, pointy, color = myCluster$clus)) + geom_point()