The KFC toy problem: perspectives from four job roles

There is an interesting question —

There are 5 different types of toys at a KFC restaurant. If you go there, you will get one toy randomly. How many times do you need to go to KFC in order to get all 5 toys?

The question is about probabilistic analysis. Different professionals, such as a business analyst, a statistical programmer, a mathematician and a software developer, will have different thinking pathway to solve this problem. Let's see what they would think.

1. Business Analyst

A business analyst will tend to do scenario analysis at the first step.

Best-case scenario:

Assume I am so lucky that each time I visit KFC and get a different toy, then I only need 5 times to have all the five toys. The minimum number is 5.

Worst-case scenario:

To get a different toy I need to go to KFC five times. If there is not the toy I want, I will come back with empty hand. Thus, I will have to go to KFC 5*5 = 25 times. Of course, this scenario never happens.

OK. Then the mean times I need to go to KFC seems to be in a range [5, 25). Let's give the simplest try — use the (5+25)/2 to get 15 times. The number is not accurate but at least we have an estimate.

2. Statistical Programmer

As a brute-force tool, simulation is the instant thought for a statistical programmer. Let the computer randomly create 10,000 trials — say a person plays the game 10,000 times. After averaging the results, the computer eventually tells the expected times to get the 5 toys.

I modified the SAS code from a post by Rick Wicklin, and set the maximum number of visits to KFC as 32. After 10,000 runs, the mean is 11.37. The hunch tells that this number should be quite close.

************(1)Simulate 10000 trials**************************;
proc iml;
   K = 5; /* number of toys */
   L = 32; /* max visits per trial */
   /* generate NSim trials of L visits */
   NSim = 10000;
    x = j(Nsim, L);
   call randseed(12345);
   call randgen(x, "Uniform");
   x = ceil(K*x); /* integers in [1,K] */
   /* record the visit when 5 toys are taken */
   c = j(NSim,1,0); /** allocate */
    do i = 1 to NSim;
      do j = 5 to L;
         rowUnique = countunique(x[i, 1:j]);
         if rowUnique = 5 then do;
            c[i, 1] = j;
            goto skip;
         end;
      end;
      skip:
    end;
   /* output the result */
   create _1 from c;
      append from c;
   close _1;
;quit;

data _2;
   set _1;
   /* remove the trials that didn't get 5 toys in 32 visits */
   where col1 >= 5;
 run;
************(2)Show the results*******************************;
proc sgplot;
   density col1;
   density col1 / type = kernel;
run;

proc means;
run;

3. Mathematician

Actually this question is a variant of Coupon collector's problem. The mean/expectation and standard deviation can be derived directly by the formulas. The final expectation should be 5*(1+1/2+1/3+1/4+1/5) = 11.41. This is the answer.

4. Software Developer

A software developer considers time complexity and space complexity first. When N is approaching infinity, the question is similar to a merge sorting. Given a merge sort is O(nlog(n)), the expected times must be greater than 5*ln(5) = 8.05. At least this number will be a lower bound for this question.

Use MongoDB as a JSON factory

MongoDB is a persistent data store for JSON formatted data, which seems like an ideal middleware between the data tier software and the web. With MongoDB, Javascript's Map/Reduce functionality makes many trivial jobs particularly easy, such as translate an object to an array. For example, we can produce a bubble plot with Highcharts.js and the SASHEP.IRIS dataset in SAS very quickly.

Step 1: push SAS dataset toward MongoDB

First, let's push the SASHEP.IRIS dataset from SAS to MongoDB using thesas2mongo macro.

%sas2mongo(data = sashelp.iris, dbname = demo, collname = iris, tmpfolder = c:\tmp, mongofolder =c:\mongodb\bin);

Step 2: make the JSON file

Under the MongoDB shell, we could use Javascript'smap function to transform data to the desired structure.

var species = ["Setosa", "Versicolor", "Virginica"];
for (i=0; i<3; i++) {
   var z = db.iris.find({Species: species[i]}, {SepalLength:1,SepalWidth:1, PetalWidth:1, _id:0}).toArray().map(function(d){return [d.SepalLength, d.SepalWidth, d.PetalWidth]});
   print("{name:", JSON.stringify(species[i]), ", data:", JSON.stringify(z), "},");
};

Step 3: finalize Highcharts.js

The link of the final plots is here.

The plots are demonstrated into a Bootstrap 3 framework. One great advantages for SVG is that it is responsive to the devices‘ screens, which is especially friendly to mobile. The PNG or JPG formatted images can invoke Bootstrap’s responsive library to have the same effect.

A SAS macro that exports data to MongoDB

MongoDB is possibly the most popular NoSQL data store. To bypass schema and constraint, I feel quite convenient to implement MongoDB as buffer to accompany current RDBMS .Also it is straightforward to use MongoDB and other tools (MEAN) to build some simple web apps for statistics presentation.

Neither SAS nor 10gen so far published any SAS-MongoDB driver. However, the table-like dataset in SAS can be transformed to CSV by PROC EXPORT or DATA Step. MongoDB has a nice API mongoimport that easily accepts CSV formatted data. I write a macro in SAS below to transport data from SAS to MongoDB. The speed is quite fast.

****************(1) MODULE-BUILDING STEP********************************;
%macro sas2mongo(data =, dbname = , collname =, tmpfolder =, mongofolder = );
   /*************************************************************************
   *  MACRO:      sas2mongo()
   *  GOAL:       output a dataset in SAS to a collection in MongoDB
   *  PARAMETERS: data       = SAS dataset to export
   *              dbname     = database name in MongoDB
   *              collname   = collection name in MongoDB
   *              tmpfolder  = Windows directory for temporary file exchange
   *              mongofolder= bin directory where MongoDB was installed
   *************************************************************************/
    proc export data=&data outfile="&tmpfolder.\tmp.csv" dbms=csv replace;
    run;
    options noxsync noxwait;
    %put the execuated command is: &mongofolder\mongoimport.exe -d 
        &dbname -c &collname --type csv --file &tmpfolder.\tmp.csv --headerline;
    x "&mongofolder\mongoimport.exe -d &dbname -c &collname --type 
        csv --file &tmpfolder.\tmp.csv --headerline";
%mend;

****************(2) TESTING STEP****************************************;
%sas2mongo(data = sashelp.class, dbname = demo, collname = class,
    tmpfolder = c:\tmp, mongofolder =c:\mongodb\bin);

Then I run commands in Mongo shell. It works just well.

use demo;
db.class.find();

Bubble plot by SAS and Highcharts.js

Bubble plot is a nice data visualization choice for three dimensional numeric variables. It seems quite popular on web and documents.

Static plotting by SAS

Since SAS 9.3, PROC SGPLOT provides a bubble statement, which makes a bubble plot easy. For example, the dataset SASHELP.CLASS can be quickly projected onto a bubble plot.

proc sgplot data = sashelp.class;
   title 'bubble plot by sashelp.class';
   bubble x = weight y = height size = age / group = sex transparency = 0.5;
   yaxis grid;
run;

Dynamic plotting by SAS and Highcharts.js

For the show-off on web, an interactive bubble plot above will be much more attractive. First we need to use SAS to transform the SASHELP.CLASS dataset to a nested JSON array. The link of the final dynamic plot is here.

data one;
   set sashelp.class;
   length data $20.;
   data = cats('[', weight, ',', height, ',', age, '],');
run;

proc sort data = one;
   by sex;
run;
proc transpose data = one out = two;
   by sex;
   var data;
run; 

data JSON;
   set two;
   length _tmp dataline $300.;
   _tmp = cats( of col:);
   substr(_tmp, length(_tmp), 1) = ' ';
   dataline = cats('{data:[', _tmp, '],', 'name:"', sex, '"}');
   keep dataline;
run;

One good thing is that the JSON data can be fully embedded in Highcharts.js, which doesn't require an HTTP server like D3.js. We only need to insert the data from SAS's DATA Step into Highcharts's bubble plot API. It also provides rich options for better visualization effects and convenient downloading.

$(function () {
    $('#container').highcharts({
        chart: {
            type: 'bubble',
            zoomType: 'xy'
        },
        credits: {
            text: "Demo",
            href: 'http://www.sasanalysis.com'
        },
        title: {
            text: 'Bubble plot by sashelp.class'
        },
        series: [{
            data: [
                [84, 56.5, 13],
                [98, 65.3, 13],
                [102.5, 62.8, 14],
                [84.5, 59.8, 12],
                [112.5, 62.5, 15],
                [50.5, 51.3, 11],
                [90, 64.3, 14],
                [77, 56.3, 12],
                [112, 66.5, 15]
            ],
            name: "F"
        }, {
            data: [
                [112.5, 69, 14],
                [102.5, 63.5, 14],
                [83, 57.3, 12],
                [84, 62.5, 13],
                [99.5, 59, 12],
                [150, 72, 16],
                [128, 64.8, 12],
                [133, 67, 15],
                [85, 57.5, 11],
                [112, 66.5, 15]
            ],
            name: "M"
        }]
    });
});

Conclusion

1. Besides its statistical feature, SAS is also a flexible scripting language such as creating JSON;
2. Highcharts.js is a view tier tool halfway between D3.js (open source; minimum documentation) and tableau(propriety software; company support), which allows integration with SAS or other data tier tools.

More SQL taste in SAS 9.4

Compared with SAS 9.3, the latest SAS 9.4 introduced a few new procedures for the BASE and STAT components: 7 new procedures for BASE 9.4 and 4 for STAT 12.3. 6 high-performance procedures (thanks to Dr. Wicklin's correction).

New in BASE 9.4	New in STAT 12.1	New in STAT 12.3
DELETE	ADAPTIVEREG	HPGENSELECT
DS2	QUANTLIFE	HPLOGISTIC
JSON	QUANTSELECT	HPLMIXED
PRESENV	STDRATE	HPNLMOD
STREAM		HPREG
FEDSQL		HPSPLIT
AUTHLIB

DS2 is a new SAS proprietary programming language that is appropriate for advanced data manipulation.It is exciting to see the emergence of DS2 and FEDSQL. According to SAS 9.4 DS2 Language Reference,

DS2 is a SAS programming language that is appropriate for advanced data manipulation

Contrary to the thought I had last year, DS2 or PROC DS2 is not a complied language. It seems more like a wrapper of PROC FEDSQL, which combines the capacity of SQL and the original DATA Step together. Therefore, DS2 includes many SQL's features such as subquery.

data class;
   set sashelp.class;
run;

proc datasets nolist;
   delete _:;
quit;

proc ds2 stimer;
   data _test1;
      dcl varchar(6) gender;
      method run();
         set {select name, sex from class where age > 12};
         if sex = 'M' then gender = 'Male';
         else gender = 'Female';
      end;
   enddata;
   run;
quit;

The functionality is equivalent to the SQL syntax in SAS below.

proc sql stimer;
   create table _test2 as
   select *, case when sex = 'M' then 'Male'
         else 'Female'
      end as gender
   from (select name, sex from class where age > 12)
;quit;

Additionally, DS2 supports the concept of transaction in SQL. The run statement in DS2 is equal to the COMMIT statement in SQL, while run cancel statement is comparable to SQL'sROLLBACK statement.

In conclusion, with DS2, SAS is leaning toward RDBMS in how to understand and deal with data.

Regularization adjustment for PROC SVM

SVM is a popular statistical learning method for either classification or regression. For classification, a linear classifier or a hyperplane, such as with w as weight vector and b as the bias, would label data into various categories. The geometric margin is defined as . For SVM, the maximum margin approach is equivalent to minimize . However, with the introduction of regularization to inhibit complexity, the optimization has to upgrade to minimize , where C is the regularization parameter. Eventually the solution for SVM turns out to be a quadratic optimization problem over w and ξ with the constraints of .

Since the sashelp.class dataset is extremely simple, I attempt to use its variables age and weight to predict sex, which is just for demonstration purpose. According to the plot, the data points are linearly non-separable. The kernel methods have to be applied to map the input data to a high-dimensional space so that they are linearly separable. To harness SVM in SAS, three procedures are commonly used under the license of SAS EMiner. For example, PROC DMDB is used to recode the categorical data and set up the working catalog, PROC SVM is used to build the model, and PROC SVMSCORE is applied to implement the model.

proc sgplot data=sashelp.class;
   scatter x = weight y = age / group = sex;
run;

proc dmdb batch data=sashelp.class dmdbcat=_cat out=_class;
   var weight age;
   class sex;
run;

Hard margin

If we let C be infinitely large, then all constraints will be executed. Therefore, the margin is narrowed down.

proc svm data=_class dmdbcat=_cat c=1e11 kernel=linear out= _1;
   title 'hard margin';
   ods output restab = restab1;
   var weight age ;
   target sex;
run;

The accuracy is 63.16%. Overall, the result is below.

Name	Value
Regularization Parameter C	100000000000
Classification Error (Training)	7.000000
Geometric Margin	1.624447E-10
Number of Support Vectors	17
Estimated VC Dim of Classifier	3.4494098E24
Number of Kernel Calls	74

Soft margin

On the contrary, the small C allows constraints to be easily ignored, which leads to the desired large margin.

proc svm data=_class dmdbcat=_cat kernel=linear out= _2;
   title 'soft margin';
   ods output restab = restab2;
   var weight age;
   target sex;
run;

The accuracy or miscalculation rate keeps the same, since the data is so small. In PROC SVM, without the specification, the C value is solved to be almost near zero, and the margin are huge.

Name	Value
Regularization Parameter C	0.000098161
Classification Error (Training)	7.000000
Geometric Margin	158.553426
Number of Support Vectors	18
Estimated VC Dim of Classifier	3.850370
Number of Kernel Calls	76

Conclusion

1. For the SVM procedure, except the training data, adding a validation data for the testdata option at the PROC statement could effectivley increase the C parameter and decrease the possibility of overfitting.
2. There are a few advantages for SVM over other data mining methods. First SVM is suitable for high dimension data, and more importantly the complexity can be easily controlled by the adjustment of the regularization parameter C.

Cluster analysis on a pivot table

The link of the pivot table is here

The increasing supremacy of JavaScript on both server side and client side seems a good news for those statistical workers who deal with data and model, and therefore always live in the darkness. They could eventually find a relatively easier way to show off their hard work on Web, the final destination of data. Here I show how to display the result of a cluster analysis on a web-based pivot table.

Back-end: cluster analysis

SAS has a FASTCLUS procedure, which implements a nearest centroid sorting algorithm and is similar to k-means. It has some time and space advantages over other more complicated clustering algorithms in SAS.

I still use the SASHELP.CLASS dataset and cluster the rows by weight and height. I specify 2 clusters and easily obtain the distances to the centroids by PROC FASTCLUS. The plot demonstrates thatweight=100 looks like the boundary to separate the two clusters. Next in DATA Step, I translate the SAS dataset to JSON format so that the browser can understand it.

************(1) Cluster the dataset*******;
proc fastclus data = sashelp.class maxclusters = 2 out = class;
   var height weight;
run;

proc sgplot data = class;
   scatter x = height y = weight /group = cluster;
   yaxis grid;
run;

************(2) Transform to JSON*********;
data toJSON;
   set class;
   length line $200.;
   array a[5] _numeric_; 
   array _a[5] $20.;
   do i = 1 to 5;
      _a[i] = cat('"',vname(a[i]),'":', a[i], ',');
   end;
   array b[2] name sex;
   array _b[2] $20.;
   do j = 1 to 2;
      _b[j] =cat('"',vname(b[j]),'":"', b[j], '",');
   end;
   line = cats('{', cats(of _:), '},');
   substr(line, length(line)-2, 1) = ' ';
   keep line;
run;

Front-end: pivot table

Pivot table is a nice way to present data, especially raw data. There are a few approaches to realize pivot table on web, such as Google's fusion table. Nicolas Kruchten developed a framework called PivotTable.js on github, which is very popular.

I embed the JSON data with the PivotTable.js to make the HTML file static, since the Blogger doesn't provide the function of HTTP server. The file content will be like:

Eventually we can view the cluster result on a pivot table. The audience can now interactively play with data.