Big Data Analytics 简明教程

# Install the following packages by running the following code in R.
pkgs = c('data.table', 'ggplot2', 'nycflights13', 'reshape2')
install.packages(pkgs)

library(nycflights13)
library(ggplot2)
library(data.table)
library(reshape2)

# Convert the flights data.frame to a data.table object and call it DT
DT <- as.data.table(flights)

# The data has 336776 rows and 16 columns
dim(DT)

# Take a look at the first rows
head(DT)

#   year    month  day   dep_time  dep_delay  arr_time  arr_delay  carrier
# 1: 2013     1     1      517       2         830         11       UA
# 2: 2013     1     1      533       4         850         20       UA
# 3: 2013     1     1      542       2         923         33       AA
# 4: 2013     1     1      544      -1        1004        -18       B6
# 5: 2013     1     1      554      -6         812        -25       DL
# 6: 2013     1     1      554      -4         740         12       UA

#     tailnum  flight  origin   dest    air_time   distance    hour   minute
# 1:  N14228   1545     EWR      IAH      227        1400       5       17
# 2:  N24211   1714     LGA      IAH      227        1416       5       33
# 3:  N619AA   1141     JFK      MIA      160        1089       5       42
# 4:  N804JB    725     JFK      BQN      183        1576       5       44
# 5:  N668DN    461     LGA      ATL      116        762        5       54
# 6:  N39463   1696     EWR      ORD      150        719        5       54

### Data Summarization
# Compute the mean arrival delay
DT[, list(mean_arrival_delay = mean(arr_delay, na.rm = TRUE))]
#        mean_arrival_delay
# 1:           6.895377
# Now, we compute the same value but for each carrier
mean1 = DT[, list(mean_arrival_delay = mean(arr_delay, na.rm = TRUE)),
   by = carrier]
print(mean1)
#      carrier    mean_arrival_delay
# 1:      UA          3.5580111
# 2:      AA          0.3642909
# 3:      B6          9.4579733
# 4:      DL          1.6443409
# 5:      EV         15.7964311
# 6:      MQ         10.7747334
# 7:      US          2.1295951
# 8:      WN          9.6491199
# 9:      VX          1.7644644
# 10:     FL         20.1159055
# 11:     AS         -9.9308886
# 12:     9E          7.3796692
# 13:     F9         21.9207048
# 14:     HA         -6.9152047
# 15:     YV         15.5569853
# 16:     OO         11.9310345

# Now let’s compute to means in the same line of code
mean2 = DT[, list(mean_departure_delay = mean(dep_delay, na.rm = TRUE),
   mean_arrival_delay = mean(arr_delay, na.rm = TRUE)),
   by = carrier]
print(mean2)

#       carrier    mean_departure_delay   mean_arrival_delay
# 1:      UA            12.106073          3.5580111
# 2:      AA             8.586016          0.3642909
# 3:      B6            13.022522          9.4579733
# 4:      DL             9.264505          1.6443409
# 5:      EV            19.955390         15.7964311
# 6:      MQ            10.552041         10.7747334
# 7:      US             3.782418          2.1295951
# 8:      WN            17.711744          9.6491199
# 9:      VX            12.869421          1.7644644
# 10:     FL            18.726075         20.1159055
# 11:     AS             5.804775         -9.9308886
# 12:     9E            16.725769          7.3796692
# 13:     F9            20.215543         21.9207048
# 14:     HA             4.900585         -6.9152047
# 15:     YV            18.996330         15.5569853
# 16:     OO            12.586207         11.9310345

### Create a new variable called gain
# this is the difference between arrival delay and departure delay
DT[, gain:= arr_delay - dep_delay]

# Compute the median gain per carrier
median_gain = DT[, median(gain, na.rm = TRUE), by = carrier]
print(median_gain)

library(nycflights13)
library(ggplot2)
library(data.table)
library(reshape2)

# Using the code from the previous section
# This computes the mean arrival and departure delays by carrier.
DT <- as.data.table(flights)
mean2 = DT[, list(mean_departure_delay = mean(dep_delay, na.rm = TRUE),
   mean_arrival_delay = mean(arr_delay, na.rm = TRUE)),
   by = carrier]

# In order to plot data in R usign ggplot, it is normally needed to reshape the data
# We want to have the data in long format for plotting with ggplot
dt = melt(mean2, id.vars = ’carrier’)

# Take a look at the first rows
print(head(dt))

# Take a look at the help for ?geom_point and geom_line to find similar examples
# Here we take the carrier code as the x axis
# the value from the dt data.table goes in the y axis

# The variable column represents the color
p = ggplot(dt, aes(x = carrier, y = value, color = variable, group = variable)) +
   geom_point() + # Plots points
   geom_line() + # Plots lines
   theme_bw() + # Uses a white background
   labs(list(title = 'Mean arrival and departure delay by carrier',
      x = 'Carrier', y = 'Mean delay'))
print(p)

# Save the plot to disk
ggsave('mean_delay_by_carrier.png', p,
   width = 10.4, height = 5.07)

# Install the package corrplot by running
install.packages('corrplot')

# then load the library
library(corrplot)

# Load the following libraries
library(nycflights13)
library(ggplot2)
library(data.table)
library(reshape2)

# We will continue working with the flights data
DT <- as.data.table(flights)
head(DT) # take a look

# We select the numeric variables after inspecting the first rows.
numeric_variables = c('dep_time', 'dep_delay',
   'arr_time', 'arr_delay', 'air_time', 'distance')

# Select numeric variables from the DT data.table
dt_num = DT[, numeric_variables, with = FALSE]

# Compute the correlation matrix of dt_num
cor_mat = cor(dt_num, use = "complete.obs")

print(cor_mat)
### Here is the correlation matrix
#              dep_time   dep_delay   arr_time   arr_delay    air_time    distance
# dep_time   1.00000000  0.25961272 0.66250900  0.23230573 -0.01461948 -0.01413373
# dep_delay  0.25961272  1.00000000 0.02942101  0.91480276 -0.02240508 -0.02168090
# arr_time   0.66250900  0.02942101 1.00000000  0.02448214  0.05429603  0.04718917
# arr_delay  0.23230573  0.91480276 0.02448214  1.00000000 -0.03529709 -0.06186776
# air_time  -0.01461948 -0.02240508 0.05429603 -0.03529709  1.00000000  0.99064965
# distance  -0.01413373 -0.02168090 0.04718917 -0.06186776  0.99064965  1.00000000

# We can display it visually to get a better understanding of the data
corrplot.mixed(cor_mat, lower = "circle", upper = "ellipse")

# save it to disk
png('corrplot.png')
print(corrplot.mixed(cor_mat, lower = "circle", upper = "ellipse"))
dev.off()

source('data_visualization.R')
### Analyzing Distributions using box-plots
# The following shows the distance as a function of the carrier

p = ggplot(DT, aes(x = carrier, y = distance, fill = carrier)) + # Define the carrier
   in the x axis and distance in the y axis
   geom_box-plot() + # Use the box-plot geom
   theme_bw() + # Leave a white background - More in line with tufte's
      principles than the default
   guides(fill = FALSE) + # Remove legend
   labs(list(title = 'Distance as a function of carrier', # Add labels
      x = 'Carrier', y = 'Distance'))
p
# Save to disk
png(‘boxplot_carrier.png’)
print(p)
dev.off()

# Let's add now another variable, the month of each flight
# We will be using facet_wrap for this
p = ggplot(DT, aes(carrier, distance, fill = carrier)) +
   geom_box-plot() +
   theme_bw() +
   guides(fill = FALSE) +
   facet_wrap(~month) + # This creates the trellis plot with the by month variable
   labs(list(title = 'Distance as a function of carrier by month',
      x = 'Carrier', y = 'Distance'))
p
# The plot shows there aren't clear differences between distance in different months

# Save to disk
png('boxplot_carrier_by_month.png')
print(p)
dev.off()

# Create a vector of numbers
numbers = c(1, 2, 3, 4, 5)
print(numbers)

# [1] 1 2 3 4 5
# Create a vector of letters
ltrs = c('a', 'b', 'c', 'd', 'e')
# [1] "a" "b" "c" "d" "e"

# Concatenate both
mixed_vec = c(numbers, ltrs)
print(mixed_vec)
# [1] "1" "2" "3" "4" "5" "a" "b" "c" "d" "e"

### Evaluate the data types using class

### One dimensional objects
# Integer vector
num = 1:10
class(num)
# [1] "integer"

# Numeric vector, it has a float, 10.5
num = c(1:10, 10.5)
class(num)
# [1] "numeric"

# Character vector
ltrs = letters[1:10]
class(ltrs)
# [1] "character"

# Factor vector
fac = as.factor(ltrs)
class(fac)
# [1] "factor"

# Matrix
M = matrix(1:12, ncol = 4)
#      [,1] [,2] [,3] [,4]
# [1,]    1    4    7   10
# [2,]    2    5    8   11
# [3,]    3    6    9   12
lM = matrix(letters[1:12], ncol = 4)
#     [,1] [,2] [,3] [,4]
# [1,] "a"  "d"  "g"  "j"
# [2,] "b"  "e"  "h"  "k"
# [3,] "c"  "f"  "i"  "l"

# Coerces the numbers to character
# cbind concatenates two matrices (or vectors) in one matrix
cbind(M, lM)
#     [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8]
# [1,] "1"  "4"  "7"  "10" "a"  "d"  "g"  "j"
# [2,] "2"  "5"  "8"  "11" "b"  "e"  "h"  "k"
# [3,] "3"  "6"  "9"  "12" "c"  "f"  "i"  "l"

class(M)
# [1] "matrix"
class(lM)
# [1] "matrix"

# data.frame
# One of the main objects of R, handles different data types in the same object.
# It is possible to have numeric, character and factor vectors in the same data.frame

df = data.frame(n = 1:5, l = letters[1:5])
df
#   n l
# 1 1 a
# 2 2 b
# 3 3 c
# 4 4 d
# 5 5 e

# Let's create a data.frame
df = data.frame(numbers = 1:26, letters)
head(df)
#      numbers  letters
# 1       1       a
# 2       2       b
# 3       3       c
# 4       4       d
# 5       5       e
# 6       6       f

# str gives the structure of a data.frame, it’s a good summary to inspect an object
str(df)
#   'data.frame': 26 obs. of  2 variables:
#   $ numbers: int  1 2 3 4 5 6 7 8 9 10 ...
#   $ letters: Factor w/ 26 levels "a","b","c","d",..: 1 2 3 4 5 6 7 8 9 10 ...

# The latter shows the letters character vector was coerced as a factor.
# This can be explained by the stringsAsFactors = TRUE argumnet in data.frame
# read ?data.frame for more information

class(df)
# [1] "data.frame"

### Indexing
# Get the first row
df[1, ]
#     numbers  letters
# 1       1       a

# Used for programming normally - returns the output as a list
df[1, , drop = TRUE]
# $numbers
# [1] 1
#
# $letters
# [1] a
# Levels: a b c d e f g h i j k l m n o p q r s t u v w x y z

# Get several rows of the data.frame
df[5:7, ]
#      numbers  letters
# 5       5       e
# 6       6       f
# 7       7       g

### Add one column that mixes the numeric column with the factor column
df$mixed = paste(df$numbers, df$letters, sep = ’’)

str(df)
# 'data.frame': 26 obs. of  3 variables:
# $ numbers: int  1 2 3 4 5 6 7 8 9 10 ...
# $ letters: Factor w/ 26 levels "a","b","c","d",..: 1 2 3 4 5 6 7 8 9 10 ...
# $ mixed  : chr  "1a" "2b" "3c" "4d" ...

### Get columns
# Get the first column
df[, 1]
# It returns a one dimensional vector with that column

# Get two columns
df2 = df[, 1:2]
head(df2)

#      numbers  letters
# 1       1       a
# 2       2       b
# 3       3       c
# 4       4       d
# 5       5       e
# 6       6       f

# Get the first and third columns
df3 = df[, c(1, 3)]
df3[1:3, ]

#      numbers  mixed
# 1       1     1a
# 2       2     2b
# 3       3     3c

### Index columns from their names
names(df)
# [1] "numbers" "letters" "mixed"
# This is the best practice in programming, as many times indeces change, but
variable names don’t
# We create a variable with the names we want to subset
keep_vars = c("numbers", "mixed")
df4 = df[, keep_vars]

head(df4)
#      numbers  mixed
# 1       1     1a
# 2       2     2b
# 3       3     3c
# 4       4     4d
# 5       5     5e
# 6       6     6f

### subset rows and columns
# Keep the first five rows
df5 = df[1:5, keep_vars]
df5

#      numbers  mixed
# 1       1     1a
# 2       2     2b
# 3       3     3c
# 4       4     4d
# 5       5     5e

# subset rows using a logical condition
df6 = df[df$numbers < 10, keep_vars]
df6

#      numbers  mixed
# 1       1     1a
# 2       2     2b
# 3       3     3c
# 4       4     4d
# 5       5     5e
# 6       6     6f
# 7       7     7g
# 8       8     8h
# 9       9     9i

# Install the sqldf package
install.packages('sqldf')

# load the library
library('sqldf')
library(nycflights13)

# We will be working with the fligths dataset in order to introduce SQL

# Let’s take a look at the table
str(flights)
# Classes 'tbl_d', 'tbl' and 'data.frame': 336776 obs. of  16 variables:

# $ year     : int  2013 2013 2013 2013 2013 2013 2013 2013 2013 2013 ...
# $ month    : int  1 1 1 1 1 1 1 1 1 1 ...
# $ day      : int  1 1 1 1 1 1 1 1 1 1 ...
# $ dep_time : int  517 533 542 544 554 554 555 557 557 558 ...
# $ dep_delay: num  2 4 2 -1 -6 -4 -5 -3 -3 -2 ...
# $ arr_time : int  830 850 923 1004 812 740 913 709 838 753 ...
# $ arr_delay: num  11 20 33 -18 -25 12 19 -14 -8 8 ...
# $ carrier  : chr  "UA" "UA" "AA" "B6" ...

# $ tailnum  : chr  "N14228" "N24211" "N619AA" "N804JB" ...
# $ flight   : int  1545 1714 1141 725 461 1696 507 5708 79 301 ...
# $ origin   : chr  "EWR" "LGA" "JFK" "JFK" ...
# $ dest     : chr  "IAH" "IAH" "MIA" "BQN" ...
# $ air_time : num  227 227 160 183 116 150 158 53 140 138 ...
# $ distance : num  1400 1416 1089 1576 762 ...
# $ hour     : num  5 5 5 5 5 5 5 5 5 5 ...
# $ minute   : num  17 33 42 44 54 54 55 57 57 58 ...

### SELECT statement
ej1 = sqldf("
   SELECT
   dep_time
   ,dep_delay
   ,arr_time
   ,carrier
   ,tailnum
   FROM
   flights
")

head(ej1)
#    dep_time   dep_delay  arr_time  carrier  tailnum
# 1      517         2      830      UA       N14228
# 2      533         4      850      UA       N24211
# 3      542         2      923      AA       N619AA
# 4      544        -1     1004      B6       N804JB
# 5      554        -6      812      DL       N668DN
# 6      554        -4      740      UA       N39463

# In R we can use SQL with the sqldf function. It works exactly the same as in
a database

# The data.frame (in this case flights) represents the table we are querying
and goes in the FROM statement
# We can also compute new variables in the select statement using the syntax:

# old_variables as new_variable
ej2 = sqldf("
   SELECT
   arr_delay - dep_delay as gain,
   carrier
   FROM
   flights
")

ej2[1:5, ]
#    gain   carrier
# 1    9      UA
# 2   16      UA
# 3   31      AA
# 4  -17      B6
# 5  -19      DL

### GROUP BY

# Computing the average
ej3 = sqldf("
  SELECT
   avg(arr_delay) as mean_arr_delay,
   avg(dep_delay) as mean_dep_delay,
   carrier
   FROM
   flights
   GROUP BY
   carrier
")

#    mean_arr_delay   mean_dep_delay carrier
# 1       7.3796692      16.725769      9E
# 2       0.3642909       8.586016      AA
# 3      -9.9308886       5.804775      AS
# 4       9.4579733      13.022522      B6
# 5       1.6443409       9.264505      DL
# 6      15.7964311      19.955390      EV
# 7      21.9207048      20.215543      F9
# 8      20.1159055      18.726075      FL
# 9      -6.9152047       4.900585      HA
# 10     10.7747334      10.552041      MQ
# 11     11.9310345      12.586207      OO
# 12      3.5580111      12.106073      UA
# 13      2.1295951       3.782418      US
# 14      1.7644644      12.869421      VX
# 15      9.6491199      17.711744      WN
# 16     15.5569853      18.996330      YV

# Other aggregations
ej4 = sqldf("
   SELECT
   avg(arr_delay) as mean_arr_delay,
   min(dep_delay) as min_dep_delay,
   max(dep_delay) as max_dep_delay,
   carrier
   FROM
   flights
   GROUP BY
   carrier
")

# We can compute the minimun, mean, and maximum values of a numeric value
ej4
#      mean_arr_delay    min_dep_delay   max_dep_delay   carrier
# 1       7.3796692           -24           747          9E
# 2       0.3642909           -24          1014          AA
# 3      -9.9308886           -21           225          AS
# 4       9.4579733           -43           502          B6
# 5       1.6443409           -33           960         DL
# 6      15.7964311           -32           548         EV
# 7      21.9207048           -27           853         F9
# 8      20.1159055           -22           602         FL
# 9      -6.9152047           -16          1301         HA
# 10     10.7747334           -26          1137         MQ
# 11     11.9310345           -14           154         OO
# 12      3.5580111           -20           483         UA
# 13      2.1295951           -19           500         US
# 14      1.7644644           -20           653         VX
# 15      9.6491199           -13           471         WN
# 16     15.5569853           -16           387         YV

### We could be also interested in knowing how many observations each carrier has
ej5 = sqldf("
   SELECT
   carrier, count(*) as count
   FROM
   flights
   GROUP BY
   carrier
")

ej5
#      carrier  count
# 1       9E    18460
# 2       AA   32729
# 3       AS   714
# 4       B6   54635
# 5       DL   48110
# 6       EV   54173
# 7       F9   685
# 8       FL   3260
# 9       HA   342
# 10      MQ   26397
# 11      OO   32
# 12      UA   58665
# 13      US   20536
# 14      VX   5162
# 15      WN   12275
# 16      YV   601

# Let’s create two tables: A and B to demonstrate joins.
A = data.frame(c1 = 1:4, c2 = letters[1:4])
B = data.frame(c1 = c(2,4,5,6), c2 = letters[c(2:5)])

A
# c1 c2
# 1  a
# 2  b
# 3  c
# 4  d

B
# c1 c2
# 2  b
# 4  c
# 5  d
# 6  e

### INNER JOIN
# This means to match the observations of the column we would join the tables by.
inner = sqldf("
   SELECT
   A.c1, B.c2
   FROM
   A INNER JOIN B
   ON A.c1 = B.c1
")

# Only the rows that match c1 in both A and B are returned
inner
# c1 c2
#  2  b
#  4  c

### LEFT OUTER JOIN
# the left outer join, sometimes just called left join will return the
# first all the values of the column used from the A table
left = sqldf("
  SELECT
   A.c1, B.c2
  FROM
   A LEFT OUTER JOIN B
   ON A.c1 = B.c1
")

# Only the rows that match c1 in both A and B are returned
left
#   c1    c2
#    1  <NA>
#    2    b
#    3  <NA>
#    4    c

# Install these packages
pkgs = c("klaR", "caret", "ElemStatLearn")
install.packages(pkgs)
library('ElemStatLearn')
library("klaR")
library("caret")

# Split the data in training and testing
inx = sample(nrow(spam), round(nrow(spam) * 0.9))
train = spam[inx,]
test = spam[-inx,]

# Define a matrix with features, X_train
# And a vector with class labels, y_train
X_train = train[,-58]
y_train = train$spam
X_test = test[,-58]
y_test = test$spam
# Train the model
nb_model = train(X_train, y_train, method = 'nb',
   trControl = trainControl(method = 'cv', number = 3))

# Compute
preds = predict(nb_model$finalModel, X_test)$class
tbl = table(y_test, yhat = preds)
sum(diag(tbl)) / sum(tbl)
# 0.7217391

library(ggplot2)
# Prepare Data
data = mtcars

# We need to scale the data to have zero mean and unit variance
data <- scale(data)

# Determine number of clusters
wss <- (nrow(data)-1)*sum(apply(data,2,var))
for (i in 2:dim(data)[2]) {
   wss[i] <- sum(kmeans(data, centers = i)$withinss)
}

# Plot the clusters
plot(1:dim(data)[2], wss, type = "b", xlab = "Number of Clusters",
   ylab = "Within groups sum of squares")

# K-Means Cluster Analysis
fit <- kmeans(data, 5) # 5 cluster solution

# get cluster means
aggregate(data,by = list(fit$cluster),FUN = mean)

# append cluster assignment
data <- data.frame(data, fit$cluster)

# Load the library for doing association rules
# install.packages(’arules’)
library(arules)

# Data preprocessing
data("AdultUCI")
AdultUCI[1:2,]
AdultUCI[["fnlwgt"]] <- NULL
AdultUCI[["education-num"]] <- NULL

AdultUCI[[ "age"]] <- ordered(cut(AdultUCI[[ "age"]], c(15,25,45,65,100)),
   labels = c("Young", "Middle-aged", "Senior", "Old"))
AdultUCI[[ "hours-per-week"]] <- ordered(cut(AdultUCI[[ "hours-per-week"]],
   c(0,25,40,60,168)), labels = c("Part-time", "Full-time", "Over-time", "Workaholic"))
AdultUCI[[ "capital-gain"]] <- ordered(cut(AdultUCI[[ "capital-gain"]],
   c(-Inf,0,median(AdultUCI[[ "capital-gain"]][AdultUCI[[ "capitalgain"]]>0]),Inf)),
   labels = c("None", "Low", "High"))
AdultUCI[[ "capital-loss"]] <- ordered(cut(AdultUCI[[ "capital-loss"]],
   c(-Inf,0, median(AdultUCI[[ "capital-loss"]][AdultUCI[[ "capitalloss"]]>0]),Inf)),
   labels = c("none", "low", "high"))

# Convert the data into a transactions format
Adult <- as(AdultUCI, "transactions")
Adult
# transactions in sparse format with
# 48842 transactions (rows) and
# 115 items (columns)

summary(Adult)
# Plot frequent item-sets
itemFrequencyPlot(Adult, support = 0.1, cex.names = 0.8)

# generate rules
min_support = 0.01
confidence = 0.6
rules <- apriori(Adult, parameter = list(support = min_support, confidence = confidence))

rules
inspect(rules[100:110, ])
# lhs                             rhs                      support     confidence  lift
# {occupation = Farming-fishing} => {sex = Male}        0.02856148  0.9362416   1.4005486
# {occupation = Farming-fishing} => {race = White}      0.02831579  0.9281879   1.0855456
# {occupation = Farming-fishing} => {native-country     0.02671881  0.8758389   0.9759474
                                       = United-States}

# Install the party package
# install.packages('party')
library(party)
library(ggplot2)

head(diamonds)
# We will predict the cut of diamonds using the features available in the
diamonds dataset.
ct = ctree(cut ~ ., data = diamonds)

# plot(ct, main="Conditional Inference Tree")
# Example output
# Response:  cut
# Inputs:  carat, color, clarity, depth, table, price, x, y, z

# Number of observations:  53940
#
# 1) table <= 57; criterion = 1, statistic = 10131.878
#   2) depth <= 63; criterion = 1, statistic = 8377.279
#     3) table <= 56.4; criterion = 1, statistic = 226.423
#       4) z <= 2.64; criterion = 1, statistic = 70.393
#         5) clarity <= VS1; criterion = 0.989, statistic = 10.48
#           6) color <= E; criterion = 0.997, statistic = 12.829
#             7)*  weights = 82
#           6) color > E

#Table of prediction errors
table(predict(ct), diamonds$cut)
#            Fair  Good Very Good Premium Ideal
# Fair       1388   171        17       0    14
# Good        102  2912       499      26    27
# Very Good    54   998      3334     249   355
# Premium      44   711      5054   11915  1167
# Ideal        22   114      3178    1601 19988
# Estimated class probabilities
probs = predict(ct, newdata = diamonds, type = "prob")
probs = do.call(rbind, probs)
head(probs)

library(ElemStatLearn)
head(spam)

# Split dataset in training and testing
inx = sample(nrow(spam), round(nrow(spam) * 0.8))
train = spam[inx,]
test = spam[-inx,]

# Fit regression model
fit = glm(spam ~ ., data = train, family = binomial())
summary(fit)

# Call:
#   glm(formula = spam ~ ., family = binomial(), data = train)
#

# Deviance Residuals:
#   Min       1Q   Median       3Q      Max
# -4.5172  -0.2039   0.0000   0.1111   5.4944
# Coefficients:
# Estimate Std. Error z value Pr(>|z|)
# (Intercept) -1.511e+00  1.546e-01  -9.772  < 2e-16 ***
# A.1         -4.546e-01  2.560e-01  -1.776 0.075720 .
# A.2         -1.630e-01  7.731e-02  -2.108 0.035043 *
# A.3          1.487e-01  1.261e-01   1.179 0.238591
# A.4          2.055e+00  1.467e+00   1.401 0.161153
# A.5          6.165e-01  1.191e-01   5.177 2.25e-07 ***
# A.6          7.156e-01  2.768e-01   2.585 0.009747 **
# A.7          2.606e+00  3.917e-01   6.652 2.88e-11 ***
# A.8          6.750e-01  2.284e-01   2.955 0.003127 **
# A.9          1.197e+00  3.362e-01   3.559 0.000373 ***
# Signif. codes:  0 *** 0.001 ** 0.01 * 0.05 . 0.1  1

### Make predictions
preds = predict(fit, test, type = ’response’)
preds = ifelse(preds > 0.5, 1, 0)
tbl = table(target = test$spam, preds)
tbl

#         preds
# target    0   1
# email   535  23
# spam     46 316
sum(diag(tbl)) / sum(tbl)
# 0.925

library(tm)
library(data.table)

source('text_analytics/text_analytics_functions.R')
data = fread('text_analytics/data/profiles.txt')
rate = as.numeric(data$rate)
keep = !is.na(rate)
rate = rate[keep]

### Make bag of words of title and body
X_all = bag_words(data$user_skills[keep])
X_all = removeSparseTerms(X_all, 0.999)
X_all

# <<DocumentTermMatrix (documents: 389, terms: 1422)>>
#   Non-/sparse entries: 4057/549101
# Sparsity           : 99%
# Maximal term length: 80
# Weighting          : term frequency - inverse document frequency (normalized) (tf-idf)

### Make a sparse matrix with all the data
X_all <- as_sparseMatrix(X_all)

train_inx = 1:200
X_train = X_all[train_inx, ]
y_train = rate[train_inx]
X_test = X_all[-train_inx, ]
y_test = rate[-train_inx]

# Train a regression model
library(glmnet)
fit <- cv.glmnet(x = X_train, y = y_train,
   family = 'gaussian', alpha = 1,
   nfolds = 3, type.measure = 'mae')
plot(fit)

# Make predictions
predictions = predict(fit, newx = X_test)
predictions = as.vector(predictions[,1])
head(predictions)

# 36.23598 36.43046 51.69786 26.06811 35.13185 37.66367
# We can compute the mean absolute error for the test data
mean(abs(y_test - predictions))
# 15.02175

python csv_to_vowpal_wabbit.py

vw train_titanic.vw -f model.vw --binary --passes 20 -c -q ff --sgd --l1
0.00000001 --l2 0.0000001 --learning_rate 0.5 --loss_function logistic

-loss_function logistic
creating quadratic features for pairs: ff
using l1 regularization = 1e-08
using l2 regularization = 1e-07

final_regressor = model.vw
Num weight bits = 18
learning rate = 0.5
initial_t = 1
power_t = 0.5
decay_learning_rate = 1
using cache_file = train_titanic.vw.cache
ignoring text input in favor of cache input
num sources = 1

average    since         example   example  current  current  current
loss       last          counter   weight    label   predict  features
0.000000   0.000000          1      1.0    -1.0000   -1.0000       57
0.500000   1.000000          2      2.0     1.0000   -1.0000       57
0.250000   0.000000          4      4.0     1.0000    1.0000       57
0.375000   0.500000          8      8.0    -1.0000   -1.0000       73
0.625000   0.875000         16     16.0    -1.0000    1.0000       73
0.468750   0.312500         32     32.0    -1.0000   -1.0000       57
0.468750   0.468750         64     64.0    -1.0000    1.0000       43
0.375000   0.281250        128    128.0     1.0000   -1.0000       43
0.351562   0.328125        256    256.0     1.0000   -1.0000       43
0.359375   0.367188        512    512.0    -1.0000    1.0000       57
0.274336   0.274336       1024   1024.0    -1.0000   -1.0000       57 h
0.281938   0.289474       2048   2048.0    -1.0000   -1.0000       43 h
0.246696   0.211454       4096   4096.0    -1.0000   -1.0000       43 h
0.218922   0.191209       8192   8192.0     1.0000    1.0000       43 h

finished run
number of examples per pass = 802
passes used = 11
weighted example sum = 8822
weighted label sum = -2288
average loss = 0.179775 h
best constant = -0.530826
best constant’s loss = 0.659128
total feature number = 427878

vw -d test_titanic.vw -t -i model.vw -p predictions.txt

# Read the predictions
preds = fread('vw/predictions.txt')

# Define the sigmoid function
sigmoid = function(x) {
   1 / (1 + exp(-x))
}
probs = sigmoid(preds[[1]])

# Generate class labels
preds = ifelse(probs > 0.5, 1, 0)
head(preds)
# [1] 0 1 0 0 1 0