Showing posts with label Machine Learning. Show all posts
Showing posts with label Machine Learning. Show all posts

Sunday, March 5, 2017

Bloom filter , Probabilistic Cache and much more

Probabilistic algorithms first make an impact as incorrect,but they are unbelievably efficient.
Also of course if probabilistic error is not on critical path you do not have to worry.
Once we had a web service here customer was asking if they send us some users.
Think they are sending 10000-20000 users on every batch,and doing so frequent.

If I hit database and check if user exists, i will query for thousand times.
If I make a cache on memory, i need to make a hashmap with item size in database.
Say I have 50K user. Cache will be 50K x user name lets say.

In Bloom Filter instead of caching every time,we compute one projected value for a hash function
and put it in a bit array we occupy much much less space and get incredible performance.

And if we use a few hash functions at one time,we get more accurate results.

This algorithm always tell exact if item exists in cache.
If it says,not in cache, there is a probability it is saying wrong.

For my scenario it is not causing any problem but maybe one unnecessary database hit.

For example I am showing my cache for Country names.
It has 227 items.

Check below address for best parameters.
http://hur.st/bloomfilter?n=227&p=0.01

For our sample parameters are below.

n = 227, p = 0.01 (1 in 100) → m = 2,176 (272B), k = 7

With only 272 Byte and 7 hash functions we can get %1 error rate.

import java.util.BitSet;
import java.util.Random;

import com.sangupta.murmur.Murmur1;
import com.sangupta.murmur.Murmur2;
import com.sangupta.murmur.Murmur3;

public class BloomFilterCache {

 public static String[] COUNTRIES = new String[] { "Afghanistan ", "Albania ", "Algeria ", "American Samoa ",
   "Andorra ", "Angola ", "Anguilla ", "Antigua & Barbuda ", "Argentina ", "Armenia ", "Aruba ", "Australia ",
   "Austria ", "Azerbaijan ", "Bahamas, The ", "Bahrain ", "Bangladesh ", "Barbados ", "Belarus ", "Belgium ",
   "Belize ", "Benin ", "Bermuda ", "Bhutan ", "Bolivia ", "Bosnia & Herzegovina ", "Botswana ", "Brazil ",
   "British Virgin Is. ", "Brunei ", "Bulgaria ", "Burkina Faso ", "Burma ", "Burundi ", "Cambodia ",
   "Cameroon ", "Canada ", "Cape Verde ", "Cayman Islands ", "Central African Rep. ", "Chad ", "Chile ",
   "China ", "Colombia ", "Comoros ", "Congo, Dem. Rep. ", "Congo, Repub. of the ", "Cook Islands ",
   "Costa Rica ", "Cote d'Ivoire ", "Croatia ", "Cuba ", "Cyprus ", "Czech Republic ", "Denmark ", "Djibouti ",
   "Dominica ", "Dominican Republic ", "East Timor ", "Ecuador ", "Egypt ", "El Salvador ",
   "Equatorial Guinea ", "Eritrea ", "Estonia ", "Ethiopia ", "Faroe Islands ", "Fiji ", "Finland ", "France ",
   "French Guiana ", "French Polynesia ", "Gabon ", "Gambia, The ", "Gaza Strip ", "Georgia ", "Germany ",
   "Ghana ", "Gibraltar ", "Greece ", "Greenland ", "Grenada ", "Guadeloupe ", "Guam ", "Guatemala ",
   "Guernsey ", "Guinea ", "Guinea-Bissau ", "Guyana ", "Haiti ", "Honduras ", "Hong Kong ", "Hungary ",
   "Iceland ", "India ", "Indonesia ", "Iran ", "Iraq ", "Ireland ", "Isle of Man ", "Israel ", "Italy ",
   "Jamaica ", "Japan ", "Jersey ", "Jordan ", "Kazakhstan ", "Kenya ", "Kiribati ", "Korea, North ",
   "Korea, South ", "Kuwait ", "Kyrgyzstan ", "Laos ", "Latvia ", "Lebanon ", "Lesotho ", "Liberia ", "Libya ",
   "Liechtenstein ", "Lithuania ", "Luxembourg ", "Macau ", "Macedonia ", "Madagascar ", "Malawi ",
   "Malaysia ", "Maldives ", "Mali ", "Malta ", "Marshall Islands ", "Martinique ", "Mauritania ",
   "Mauritius ", "Mayotte ", "Mexico ", "Micronesia, Fed. St. ", "Moldova ", "Monaco ", "Mongolia ",
   "Montserrat ", "Morocco ", "Mozambique ", "Namibia ", "Nauru ", "Nepal ", "Netherlands ",
   "Netherlands Antilles ", "New Caledonia ", "New Zealand ", "Nicaragua ", "Niger ", "Nigeria ",
   "N. Mariana Islands ", "Norway ", "Oman ", "Pakistan ", "Palau ", "Panama ", "Papua New Guinea ",
   "Paraguay ", "Peru ", "Philippines ", "Poland ", "Portugal ", "Puerto Rico ", "Qatar ", "Reunion ",
   "Romania ", "Russia ", "Rwanda ", "Saint Helena ", "Saint Kitts & Nevis ", "Saint Lucia ",
   "St Pierre & Miquelon ", "Saint Vincent and the Grenadines ", "Samoa ", "San Marino ",
   "Sao Tome & Principe ", "Saudi Arabia ", "Senegal ", "Serbia ", "Seychelles ", "Sierra Leone ",
   "Singapore ", "Slovakia ", "Slovenia ", "Solomon Islands ", "Somalia ", "South Africa ", "Spain ",
   "Sri Lanka ", "Sudan ", "Suriname ", "Swaziland ", "Sweden ", "Switzerland ", "Syria ", "Taiwan ",
   "Tajikistan ", "Tanzania ", "Thailand ", "Togo ", "Tonga ", "Trinidad & Tobago ", "Tunisia ", "Turkey ",
   "Turkmenistan ", "Turks & Caicos Is ", "Tuvalu ", "Uganda ", "Ukraine ", "United Arab Emirates ",
   "United Kingdom ", "United States ", "Uruguay ", "Uzbekistan ", "Vanuatu ", "Venezuela ", "Vietnam ",
   "Virgin Islands ", "Wallis and Futuna ", "West Bank ", "Western Sahara ", "Yemen ", "Zambia ",
   "Zimbabwe " };

 private BitSet bitset;

 int[] hashSeeds;

 int noHashFunctions;

 public boolean useSeed = true;

 public int murmurVersion = 1;

 public int getSeed(int i) {
  if (useSeed)
   return hashSeeds[i];
  else
   return i;

 }

 public BloomFilterCache(int slots, int hashFunctions) {
  bitset = new BitSet(slots);
  noHashFunctions = hashFunctions;
  Random r = new Random(System.currentTimeMillis());
  hashSeeds = new int[hashFunctions];
  for (int i = 0; i < hashFunctions; ++i) {
   hashSeeds[i] = r.nextInt();
  }

 }

 public int getHash(String value, int i) {

  if (murmurVersion == 1)
   return (int) Murmur1.hash(value.getBytes(), 4, getSeed(i));
  else if (murmurVersion == 2)
   return (int) Murmur2.hash(value.getBytes(), 4, getSeed(i));
  else
   return (int) Murmur3.hash_x86_32(value.getBytes(), 4, getSeed(i));
 }

 public void add(String value) {

  for (int i = 0; i < noHashFunctions; ++i) {
   int h = getHash(value, i);

   bitset.set(Math.abs(h) % bitset.size(), true);
  }
 }

 public boolean contains(String value) {

  for (int i = 0; i < noHashFunctions; ++i) {
   int h = getHash(value, i);

   if (!bitset.get(Math.abs(h) % bitset.size())) {
    return false;

   }
  }

  return true;
 }

 public static String generateRandomChars(String candidateChars, int length) {
  StringBuilder sb = new StringBuilder();
  Random random = new Random();
  for (int i = 0; i < length; i++) {
   sb.append(candidateChars.charAt(random.nextInt(candidateChars.length())));
  }

  return sb.toString();
 }

 public int getErrorCount(String[] testSet, String postFix, String prefix) {
  int error = 0;

  for (String s : testSet) {

   if (contains(generateRandomChars("ABCDEFGHIJKLMNOPQRSTUVWXYZ".toLowerCase(), 3) + s + "")) {
    error++;
   }
  }
  return error;
 }

 public int getErrorCountRandom(int testSize) {
  int error = 0;

  for (int i = 0; i < testSize; i++) {

   String test = generateRandomChars("ABCDEFGHIJKLMNOPQRSTUVWXYZ".toLowerCase(), 10);
   if (contains(test) ) {
    error++;
    //System.err.println(test);
   }

  }
  return error;
 }

 public static void testFor(int slots, int hashFunctions, boolean useSeed, int murmurVersion) {
  BloomFilterCache bf = new BloomFilterCache(slots, hashFunctions);
  bf.useSeed = useSeed;
  bf.murmurVersion = murmurVersion;

  for (String s : COUNTRIES) {
   bf.add(s);
  }

  System.err.println("*****Params : slots = " + slots + " no#hash = " + hashFunctions + " cardinality = "
    + bf.bitset.cardinality() + " useSeed = " + useSeed + "  murmurVersion = " + murmurVersion);
  System.err.println("Query for Japan: " + bf.contains("Japan"));
  System.err.println("Query for Dummy: " + bf.contains("Dummy"));
  System.err.println("Error Count: " + bf.getErrorCount(COUNTRIES, "", ""));
  System.err.println("Error Count Prefix: " + bf.getErrorCount(COUNTRIES, "abc", ""));
  System.err.println("Error Count Postfix:  " + bf.getErrorCount(COUNTRIES, "", "abc"));
  System.err.println("Error Count Post-Pre: " + bf.getErrorCount(COUNTRIES, "abc", "def"));
  System.err.println("Error Count Random: " + bf.getErrorCountRandom(COUNTRIES.length));

 }

 public static void main(String[] args) {


  int[] versions = new int[] { 1, 2, 3 };
  boolean[] useSeeds = new boolean[] { true, false };
  for (int v : versions) {
   System.err.println(" ---------------------- ");
   for (boolean useSeed : useSeeds) {
    testFor(10000, 7, useSeed, v);
    testFor(2224, 2, useSeed, v);
    testFor(2224, 7, useSeed, v);
    testFor(2224, 17, useSeed, v);
    testFor(582, 3, useSeed, v);
   }
  }

 }
}

Below are multiple tests with multiple parameters . You can check effects by changing parameters.
97
97
3104
96321
 ---------------------- 
*****Params : slots = 10000 no#hash = 7 cardinality = 1336 useSeed = true  murmurVersion = 1
Query for Japan: true
Query for Dummy: false
Error Count: 0
Error Count Prefix: 0
Error Count Postfix:  0
Error Count Post-Pre: 0
Error Count Random: 0
*****Params : slots = 2224 no#hash = 2 cardinality = 369 useSeed = true  murmurVersion = 1
Query for Japan: true
Query for Dummy: false
Error Count: 10
Error Count Prefix: 6
Error Count Postfix:  4
Error Count Post-Pre: 4
Error Count Random: 2
*****Params : slots = 2224 no#hash = 7 cardinality = 1067 useSeed = true  murmurVersion = 1
Query for Japan: true
Query for Dummy: false
Error Count: 3
Error Count Prefix: 3
Error Count Postfix:  1
Error Count Post-Pre: 2
Error Count Random: 1
*****Params : slots = 2224 no#hash = 17 cardinality = 1755 useSeed = true  murmurVersion = 1
Query for Japan: true
Query for Dummy: false
Error Count: 3
Error Count Prefix: 4
Error Count Postfix:  5
Error Count Post-Pre: 3
Error Count Random: 2
*****Params : slots = 582 no#hash = 3 cardinality = 402 useSeed = true  murmurVersion = 1
Query for Japan: true
Query for Dummy: false
Error Count: 61
Error Count Prefix: 57
Error Count Postfix:  53
Error Count Post-Pre: 58
Error Count Random: 54
*****Params : slots = 10000 no#hash = 7 cardinality = 1308 useSeed = false  murmurVersion = 1
Query for Japan: true
Query for Dummy: false
Error Count: 0
Error Count Prefix: 0
Error Count Postfix:  0
Error Count Post-Pre: 0
Error Count Random: 0
*****Params : slots = 2224 no#hash = 2 cardinality = 371 useSeed = false  murmurVersion = 1
Query for Japan: true
Query for Dummy: false
Error Count: 10
Error Count Prefix: 9
Error Count Postfix:  7
Error Count Post-Pre: 7
Error Count Random: 5
*****Params : slots = 2224 no#hash = 7 cardinality = 1034 useSeed = false  murmurVersion = 1
Query for Japan: true
Query for Dummy: false
Error Count: 0
Error Count Prefix: 2
Error Count Postfix:  1
Error Count Post-Pre: 2
Error Count Random: 0
*****Params : slots = 2224 no#hash = 17 cardinality = 1721 useSeed = false  murmurVersion = 1
Query for Japan: true
Query for Dummy: false
Error Count: 2
Error Count Prefix: 3
Error Count Postfix:  4
Error Count Post-Pre: 3
Error Count Random: 1
*****Params : slots = 582 no#hash = 3 cardinality = 378 useSeed = false  murmurVersion = 1
Query for Japan: true
Query for Dummy: false
Error Count: 49
Error Count Prefix: 42
Error Count Postfix:  47
Error Count Post-Pre: 33
Error Count Random: 45
 ---------------------- 
*****Params : slots = 10000 no#hash = 7 cardinality = 1339 useSeed = true  murmurVersion = 2
Query for Japan: true
Query for Dummy: false
Error Count: 0
Error Count Prefix: 0
Error Count Postfix:  0
Error Count Post-Pre: 0
Error Count Random: 0
*****Params : slots = 2224 no#hash = 2 cardinality = 382 useSeed = true  murmurVersion = 2
Query for Japan: true
Query for Dummy: false
Error Count: 9
Error Count Prefix: 6
Error Count Postfix:  6
Error Count Post-Pre: 3
Error Count Random: 4
*****Params : slots = 2224 no#hash = 7 cardinality = 1069 useSeed = true  murmurVersion = 2
Query for Japan: true
Query for Dummy: false
Error Count: 1
Error Count Prefix: 3
Error Count Postfix:  1
Error Count Post-Pre: 1
Error Count Random: 1
*****Params : slots = 2224 no#hash = 17 cardinality = 1780 useSeed = true  murmurVersion = 2
Query for Japan: true
Query for Dummy: false
Error Count: 3
Error Count Prefix: 3
Error Count Postfix:  4
Error Count Post-Pre: 3
Error Count Random: 8
*****Params : slots = 582 no#hash = 3 cardinality = 386 useSeed = true  murmurVersion = 2
Query for Japan: true
Query for Dummy: false
Error Count: 50
Error Count Prefix: 43
Error Count Postfix:  48
Error Count Post-Pre: 40
Error Count Random: 53
*****Params : slots = 10000 no#hash = 7 cardinality = 1348 useSeed = false  murmurVersion = 2
Query for Japan: true
Query for Dummy: false
Error Count: 0
Error Count Prefix: 0
Error Count Postfix:  0
Error Count Post-Pre: 0
Error Count Random: 0
*****Params : slots = 2224 no#hash = 2 cardinality = 383 useSeed = false  murmurVersion = 2
Query for Japan: true
Query for Dummy: false
Error Count: 11
Error Count Prefix: 5
Error Count Postfix:  6
Error Count Post-Pre: 8
Error Count Random: 5
*****Params : slots = 2224 no#hash = 7 cardinality = 1063 useSeed = false  murmurVersion = 2
Query for Japan: true
Query for Dummy: false
Error Count: 0
Error Count Prefix: 2
Error Count Postfix:  1
Error Count Post-Pre: 1
Error Count Random: 2
*****Params : slots = 2224 no#hash = 17 cardinality = 1749 useSeed = false  murmurVersion = 2
Query for Japan: true
Query for Dummy: false
Error Count: 1
Error Count Prefix: 3
Error Count Postfix:  4
Error Count Post-Pre: 1
Error Count Random: 3
*****Params : slots = 582 no#hash = 3 cardinality = 408 useSeed = false  murmurVersion = 2
Query for Japan: true
Query for Dummy: false
Error Count: 60
Error Count Prefix: 61
Error Count Postfix:  49
Error Count Post-Pre: 47
Error Count Random: 60
 ---------------------- 
*****Params : slots = 10000 no#hash = 7 cardinality = 1322 useSeed = true  murmurVersion = 3
Query for Japan: true
Query for Dummy: false
Error Count: 0
Error Count Prefix: 0
Error Count Postfix:  0
Error Count Post-Pre: 0
Error Count Random: 0
*****Params : slots = 2224 no#hash = 2 cardinality = 359 useSeed = true  murmurVersion = 3
Query for Japan: true
Query for Dummy: false
Error Count: 18
Error Count Prefix: 22
Error Count Postfix:  12
Error Count Post-Pre: 13
Error Count Random: 14
*****Params : slots = 2224 no#hash = 7 cardinality = 1027 useSeed = true  murmurVersion = 3
Query for Japan: true
Query for Dummy: false
Error Count: 0
Error Count Prefix: 0
Error Count Postfix:  3
Error Count Post-Pre: 8
Error Count Random: 6
*****Params : slots = 2224 no#hash = 17 cardinality = 1727 useSeed = true  murmurVersion = 3
Query for Japan: true
Query for Dummy: false
Error Count: 8
Error Count Prefix: 6
Error Count Postfix:  8
Error Count Post-Pre: 3
Error Count Random: 6
*****Params : slots = 582 no#hash = 3 cardinality = 389 useSeed = true  murmurVersion = 3
Query for Japan: true
Query for Dummy: false
Error Count: 52
Error Count Prefix: 63
Error Count Postfix:  49
Error Count Post-Pre: 49
Error Count Random: 56
*****Params : slots = 10000 no#hash = 7 cardinality = 1332 useSeed = false  murmurVersion = 3
Query for Japan: true
Query for Dummy: false
Error Count: 0
Error Count Prefix: 0
Error Count Postfix:  0
Error Count Post-Pre: 0
Error Count Random: 0
*****Params : slots = 2224 no#hash = 2 cardinality = 383 useSeed = false  murmurVersion = 3
Query for Japan: true
Query for Dummy: false
Error Count: 12
Error Count Prefix: 4
Error Count Postfix:  10
Error Count Post-Pre: 5
Error Count Random: 5
*****Params : slots = 2224 no#hash = 7 cardinality = 1060 useSeed = false  murmurVersion = 3
Query for Japan: true
Query for Dummy: false
Error Count: 2
Error Count Prefix: 0
Error Count Postfix:  3
Error Count Post-Pre: 0
Error Count Random: 1
*****Params : slots = 2224 no#hash = 17 cardinality = 1782 useSeed = false  murmurVersion = 3
Query for Japan: true
Query for Dummy: false
Error Count: 4
Error Count Prefix: 6
Error Count Postfix:  7
Error Count Post-Pre: 6
Error Count Random: 5
*****Params : slots = 582 no#hash = 3 cardinality = 392 useSeed = false  murmurVersion = 3
Query for Japan: true
Query for Dummy: false
Error Count: 59
Error Count Prefix: 59
Error Count Postfix:  57
Error Count Post-Pre: 60
Error Count Random: 57

Sunday, February 19, 2017

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.

Saturday, January 28, 2017

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 
> 

Sunday, January 15, 2017

Gradient Descent for Linear Regression


I was searching for something(I do not remember) and I saw below post.
GradientDescentExample

This was a perfect post to test some parameters on gradient descent.
I opened my databricks notebook and began to play with it. I took the functions from this page but changed a bit
because of type changes in my code.

Question : Check picture below. We have points as below, is there a formula
that identifies this spread.


I used y1 = 5 * x1 + 10 + noise formula to generate this data.
So our target values are 5 and 10. ( Or a little different because of noise)

noise = np.random.normal(-3, 6, 49)

x1 = np.linspace(0, 50, 49)
y1 = 5 * x1 + 10 + noise
points = zip(x1,y1)

fig3, ax3 = plt.subplots()
#ax3.plot(x1, y1, 'k--')

ax3.plot(x1, y1, 'ro')



display(fig3)


So lets say you made an initial guess.
y = 5 * x + 3



Now lets say this is your initial guess. We must calculate how good is y = 5 * x + 3
from sum( (guess - actual)^2 ) / len ( standard formula)

What is next step? Make a better guess. How do determine you will be making a better guess.
There must be a function which will determine how your error decreases. Gradient descent function will
help you choose better values for slope and intercept.

There are some parameters you are giving to function. Learning rate and iteration count.
Learning rate is hard to understand. I checked various learning rates to see the effect.

for learning_rate = 0.00001
As you see in the graph error function is diminishing at each run. But by time
improvement is getting smaller.



for learning_rate = 0.001
It seems it is learning faster. But be careful this is a simple example and our distribution is simple with only one minimum.
If we had a complex function who has lots of convex ,concave shapes then our high/low learning rates could skip global minimum
or stuck in local minimum. Check pictures in net for these effects. There are lots of nice pictures.


You can play with parameters below and you will obtain very different results according to your parameters.

num_tests = 10
mycoef = 1
iter_count = 500

learning_rate = 0.001
initial_b = 0 # initial y-intercept guess
initial_m = 0 # initial slope guess


figError, axError = plt.subplots()

fig2, ax2 = plt.subplots()
plt.figure(1)

ax2.plot(x1, y1, 'ro')

errorList = [];


for i in range(num_tests):
  [b, m] = gradient_descent_runner(points, initial_b, initial_m, learning_rate, (i+mycoef)* iter_count)  
  
  
  x = np.linspace(0, 50, 49)
  y = x *m + b
  ax2.plot(x, y, 'k--')
  ax2.text(max(x),max(y),i)

  error = compute_error_for_line_given_points(b, m, points)
  errorList.append( error );
  
  
  axError.plot( i ,error , 'bo')
  print "After {0} iterations b = {1}, m = {2}, error = {3}".format( (i+mycoef)* iter_count, b, m, error)

Below is a result of parameters
num_tests = 10
mycoef = 1
iter_count = 500

You can play as much as you want and see the of error, slope and intercept.
Below graph seems bad because ,function performs so good from beginning and lines overlap.



Part I took from article

from numpy import *
import numpy as np
from StringIO import StringIO
import matplotlib.pyplot as plt
import numpy as np

# y = mx + b
# m is slope, b is y-intercept
def compute_error_for_line_given_points(b, m, points):
    totalError = 0
    for i in range(0, len(points)):
        x = points[i][ 0]
        y = points[i][ 1]
        totalError += (y - (m * x + b)) ** 2
    return totalError / float(len(points))

def step_gradient(b_current, m_current, points, learningRate):
    b_gradient = 0
    m_gradient = 0
    N = float(len(points))
    for i in range(0, len(points)):
        x = points[i, 0]
        y = points[i, 1]
        b_gradient += -(2/N) * (y - ((m_current * x) + b_current))
        m_gradient += -(2/N) * x * (y - ((m_current * x) + b_current))
    new_b = b_current - (learningRate * b_gradient)
    new_m = m_current - (learningRate * m_gradient)
    return [new_b, new_m]

def gradient_descent_runner(points, starting_b, starting_m, learning_rate, num_iterations):
    b = starting_b
    m = starting_m
    for i in range(num_iterations):
        b, m = step_gradient(b, m, array(points), learning_rate)
    return [b, m]

Saturday, January 7, 2017

Trying to guess if a car is damaged or not by Logistic Regresion

We had a car data collected from website. It was an advertisement website for used cars.
Car data had below properties regarding damage.

1)If there is a big damage or any damage that insurance company knows people
say it is damaged.
2)If it is a small thing,if one can make up, or already painted that are and thinks
no one understand , he does not say it is damaged.

What we have is then
Car model
Year
Price
City
Date of publishing
Last update of advertisement
Days elapsed from publishing(if a car is sold it goes from list)
Elapsed Days for selling

Think there is no data as clue. We must generate , extract, invent our data.

a)So lets think how a damaged car owner thinks
b)What changes in advertisement over time if car has small damage.(car owner is
editing advertisement over time)



1)Number of page view

There is a mean number of average page view before car owner deletes advertisement.
Lets say a no damaged car is being sold after nearly 100 page views. If a car is advertised
as not damaged and still not sold after 100 page views it can have a problem.

2)Number of change in price

At 1st owner thinks he can sell his car with a price like non-damaged cars .After a period he makes some discounts.
Probably after some calls he realizes he has to make discount. So we can generate 2 variables from here
% discount he made from 1st price
# of discounts he made.

3)Is price lower than average with same conditions.

A sense of guilt could be determined.

4)Duration that it is on sale

Total duration car is on sale.


5)Difference of duration in days from duration average sales of same car model.

There is an average duration for every combination of cars. So elapsed day
after average duration(day or week) will increase the probability of damage.



6)Number of pictures in advertisement

Probably a damaged car owner will put no picture or 1-2 pictures.Less picture could mean more probability of damage.



Monday, December 26, 2016

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

Friday, December 23, 2016

Spark NaiveBayes and Result Interpretation

***
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.
***
In samples at internet people usually try to guess if a mail is spam or not.
Below code includes codes from spark samples and some other samples.
I tried to work with spark 2 but it was not success. It is working with 1.6.
Since it did not work i played a lot, took lots of fixes from net. So code is not neat.

Lets make it much more simpler. I will list some properties and try to guess if it is
Plane or Not.
My training is
"wing wheel engine" : 1 it is plane
"wheel airbag engine" : 0 it is not plane

Steps
1)Get training set
2)Tokenize it
3)Apply hashingtf

Result of hashingtf, it generates 2vectors of words.

0 wheel airbag engine ["wheel","airbag","engine"] {"type":0,"size":20,"indices":[3,14,18],"values":[1,1,1]}
1 wing wheel engine ["wing","wheel","engine"] {"type":0,"size":20,"indices":[3,7,14],"values":[1,1,1]}


4)Train model
Result of training
Array[org.apache.spark.mllib.regression.LabeledPoint] = Array(
(8.0,[0.0,0.0,0.0,0.0,0.0,0.0,0.0,1.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,1.0,0.0]), 
(9.0,[0.0,0.0,0.0,1.0,0.0,0.0,0.0,1.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0])
)

5)Prepare test data
(0,"wing airbag")
(1,"wing airport")
(0,"wing airport") False negative(this will be guest as plane, but it is zeppelin!!)
hashingtf generates below vectors for test data

[8,wing airbag,WrappedArray(wing, airbag),(20,[7,18],[1.0,1.0])], 
[9,wing airport,WrappedArray(wing, airport),(20,[3,7],[1.0,1.0])])

7 was wing and 3 was wheel in model vectors.

6)Apply prediction testpredictionAndLabel
(0.0,0.0) I guessed as not plane, Not plane
(1.0,1.0) I guessed as plane , Plane
(1.0,0.0) I guessed as plane , Plane

7)Dump metrics, output is as below

Confusion matrix: 
1.0 1.0 
0.0 1.0 
Precision(0.0) = 1.0 
Precision(1.0) = 0.5 
Recall(0.0) = 0.5 
Recall(1.0) = 1.0 
FPR(0.0) = 0.0 
FPR(1.0) = 0.5 
F1-Score(0.0) = 0.6666666666666666 
F1-Score(1.0) = 0.6666666666666666 
Weighted precision: 0.8333333333333333 
Weighted recall: 0.6666666666666666 
Weighted F1 score: 0.6666666666666666 
Weighted false positive rate: 0.16666666666666666 
labels: Array[Double] = Array(0.0, 1.0)

Precision (0.0) :1.0 we guest 1 zero(false) , that was correct so ratio 1 / 1 = 1
Precision (1.0) :0.5 we guest 2 one(true) , 1 was correct 1 not so ratio 1 / 2 = 0.5

Recall(0.0) :0.5 we guest 1 zero(false) , there was infact 2 zeros 1 / 2 = 0.5
Recall(1.0) :1.0 we guest 1 one(true) , there was correct 1 not so ratio 1 / 1 = 1

F1-Score(0.0) = 0.6666666666666666

F1- Score = 2 x ( precision x recall ) / precison + recall.
= 2 x ( 0.5 x 1 ) / 0.5 + 1 = 2 x 0.5 / 1.5 = 0.6

From definitions :
Precision can be seen as a measure of exactness or quality, whereas recall is a measure of completeness or quantity.
High precision means that an algorithm returned substantially more relevant results than irrelevant ones, while high recall means that an algorithm returned most of the relevant results.

What does these mean.
Think in our sample we have a bigger set and we say
there are 30 planes but only 20 of them is really (among 30)
then precision is 20 / 30 = this is how well we performed on our results.
But there are items we missed.
Think in fact there were total 50 planes.
Then recall = 20 / 50 = 0.4
It is what percent of real result we returned.


High precision Low recall : we are very good at estimation but we do not cover the whole space. It means
we choose cut-off value so high.





import org.apache.spark.ml.feature.{RegexTokenizer, Tokenizer}
import org.apache.spark.ml.feature.{HashingTF, IDF}
import org.apache.spark.mllib.classification.{NaiveBayes, NaiveBayesModel}
import org.apache.spark.mllib.util.MLUtils
import org.apache.spark.mllib.linalg.Vectors
import org.apache.spark.mllib.regression.LabeledPoint
import org.apache.spark.mllib.linalg.Vector
import org.apache.spark.mllib.evaluation.MulticlassMetrics



val trainData = sqlContext.createDataFrame(Seq((0,"wheel airbag engine"),(1,"wing wheel engine"))).toDF("category","text")
    val tokenizer = new Tokenizer().setInputCol("text").setOutputCol("words")
    val wordsData = tokenizer.transform(trainData)
    val hashTF = new HashingTF().setInputCol("words").setOutputCol("features").setNumFeatures(20)
    val featureData = hashTF.transform(wordsData) 
val subFeature = featureData.select("category","features");
val df_1 = subFeature.withColumnRenamed("category","category2")
val trainDataRdd2 = df_1.withColumn("category",df_1.col("category2").cast("double")).drop("category2")


trainDataRdd2.printSchema()
val testScoreAndLabel = trainDataRdd2.select("category","features").map{ case Row(l:Double,p:Vector) => LabeledPoint(l,p) }

    val model = NaiveBayes.train(testScoreAndLabel, lambda = 1.0, modelType = "multinomial")
   //same for the test data
    val testData = sqlContext.createDataFrame(Seq((0,"wing airbag"),(1,"wing airport"),(0,"wing airport"))).toDF("category","text")
    val testWordData = tokenizer.transform(testData)
    val testFeatureData = hashTF.transform(testWordData)
    val testDataRdd = testFeatureData.select("category","features").map {
    case Row(label: Int, features: Vector) =>
    LabeledPoint(label.toDouble, Vectors.dense(features.toArray))
    }
    val testpredictionAndLabel = testDataRdd.map(p => (model.predict(p.features), p.label))


val metrics = new MulticlassMetrics(testpredictionAndLabel)
/* output F1-measure for all labels (0 and 1, negative and positive) */
metrics.labels.foreach( l => println(metrics.fMeasure(l)))
testpredictionAndLabel.take(5)
// Confusion matrix
println("Confusion matrix:")
println(metrics.confusionMatrix)


// Precision by label
val labels = metrics.labels
labels.foreach { l =>
  println(s"Precision($l) = " + metrics.precision(l))
}

// Recall by label
labels.foreach { l =>
  println(s"Recall($l) = " + metrics.recall(l))
}

// False positive rate by label
labels.foreach { l =>
  println(s"FPR($l) = " + metrics.falsePositiveRate(l))
}

// F-measure by label
labels.foreach { l =>
  println(s"F1-Score($l) = " + metrics.fMeasure(l))
}

// Weighted stats
println(s"Weighted precision: ${metrics.weightedPrecision}")
println(s"Weighted recall: ${metrics.weightedRecall}")
println(s"Weighted F1 score: ${metrics.weightedFMeasure}")
println(s"Weighted false positive rate: ${metrics.weightedFalsePositiveRate}")
  

Thursday, December 22, 2016

Spark BinaryClassificationMetrics

After finishing a LogisticRegression we can check if result is good with BinaryClassificationMetrics.
It simply takes 2 parameters.
One is score associated with your predicition.(rawPrediction column after a Logistic Regression) for example.
And other is what you guesses.
For ROC you must have a big area near 1.

What does this mean.
Suppose you are measuring if you use heater according to weather.
(of course this is obvious, we are now doing obvious case)

Say at 10 F : do not use
20 F : do not use
...
50 F : use
..
100 F : use

You see for low scores ,u do not use, but for high ones you use.
True and false it perfectly separated so I expect a perfect ROC.

ROC is a graph showing what we gain as data for calculations we did with Logistic Regression.
for example we can have 4 data for one point. If you check below it means we only use heater once on this period.
So value 10 give 3 0 and 1 1 value. This makes learning of value 10 less efficient.
Intervals must give as much as information as possible.
Purified intervals will only output 1 value for so that information gain. is so high.

( 10.0, 0.0),
( 10.0, 0.0),
( 10.0, 0.0),
( 10.0, 1.0),


val metricData= sc.parallelize(
  
   Seq( 
     ( 10.0,  0.0),
     ( 20.0,  0.0),
     ( 30.0,  0.0),
     ( 40.0,  0.0),
     ( 50.0,  0.0),
     ( 60.0,  1.0),
     ( 70.0,  1.0),
     ( 80.0,  1.0),
     ( 90.0,  1.0),     
     ( 100.0,  1.0)
    
     )
);

val metrics = new BinaryClassificationMetrics(metricData) 
println("area under the precision-recall curve: " + metrics.areaUnderPR)
println("area under the receiver operating characteristic (ROC) curve : " + metrics.areaUnderROC)
metrics.roc().collect()



Above case was so good so metrics are below.

area under the precision-recall curve: 1.0 
area under the receiver operating characteristic (ROC) curve : 0.9999999999999999 
Array[(Double, Double)] = Array((0.0,0.0), (0.0,0.2), (0.0,0.4), (0.0,0.6), (0.0,0.8), (0.0,1.0), (0.2,1.0), (0.4,1.0), (0.6,1.0), (0.8,1.0), (1.0,1.0), (1.0,1.0))







Lets preapre a bad data where distribution is useless.
Think you are measuring your ice-tea consumption according to weather.
As in above you do need have a pattern. You do not drink at 10F but you drink 20 ...
So this is near random distribution. And random distribution gives 0.5 area under curve.
It is 45 degree line. A line like that means on every probability(score) of event
I have equal info from True or False case.


val metricData= sc.parallelize(
  
   Seq( 
     ( 10.0,  0.0),
     ( 20.0,  1.0),
     ( 30.0,  0.0),
     ( 40.0,  1.0),
     ( 50.0,  0.0),
     ( 60.0,  1.0),
     ( 70.0,  0.0),
     ( 80.0,  1.0),
     ( 90.0,  0.0),     
     ( 100.0,  1.0)
    
     )
);

val metrics = new BinaryClassificationMetrics(metricData) 
println("area under the precision-recall curve: " + metrics.areaUnderPR)
println("area under the receiver operating characteristic (ROC) curve : " + metrics.areaUnderROC)
metrics.roc().collect()







Above case was so bad so metrics are below.

area under the precision-recall curve: 0.6393650793650794 
area under the receiver operating characteristic (ROC) curve : 0.6000000000000001 metrics: 
Array[(Double, Double)] = Array((0.0,0.0), (0.0,0.2), (0.2,0.2), (0.2,0.4), (0.4,0.4), (0.4,0.6), (0.6,0.6), (0.6,0.8), (0.8,0.8), (0.8,1.0), (1.0,1.0), (1.0,1.0))