Data Loading and Initial Exploration

# Load the datasets

df_business <- read_csv('business_AZ_3.csv')
df_review <- read_csv('review_AZ_3.csv')
df_tip <- read_csv('tip_AZ_3.csv')

cat("Business dataset shape:", dim(df_business)[1], "x", dim(df_business)[2], "\n")
## Business dataset shape: 6187 x 71
cat("Review dataset shape:", dim(df_review)[1], "x", dim(df_review)[2], "\n")
## Review dataset shape: 350340 x 8
cat("Tip dataset shape:", dim(df_tip)[1], "x", dim(df_tip)[2], "\n")
## Tip dataset shape: 54219 x 4
# Check business data structure
cat("\nBusiness data sample:\n")
## 
## Business data sample:
knitr::kable(head(df_business, 10), 
             caption = "Sample of Business Data")
Sample of Business Data
business_id name address city state postal_code latitude longitude stars review_count is_open attributes.ByAppointmentOnly attributes.BusinessAcceptsCreditCards attributes.BikeParking attributes.RestaurantsPriceRange2 attributes.CoatCheck attributes.RestaurantsTakeOut attributes.RestaurantsDelivery attributes.Caters attributes.WiFi attributes.BusinessParking attributes.WheelchairAccessible attributes.HappyHour attributes.OutdoorSeating attributes.HasTV attributes.RestaurantsReservations attributes.DogsAllowed attributes.Alcohol attributes.GoodForKids attributes.RestaurantsAttire attributes.Ambience attributes.RestaurantsTableService attributes.RestaurantsGoodForGroups attributes.DriveThru attributes.NoiseLevel attributes.GoodForMeal attributes.BusinessAcceptsBitcoin attributes.Smoking attributes.Music attributes.GoodForDancing attributes.AcceptsInsurance attributes.BestNights attributes.BYOB attributes.Corkage attributes.BYOBCorkage attributes.HairSpecializesIn attributes.Open24Hours attributes.RestaurantsCounterService attributes.AgesAllowed attributes.DietaryRestrictions categories hours.Monday hours.Tuesday hours.Wednesday hours.Thursday hours.Friday hours.Saturday hours.Sunday Food_Beverage Health_Fitness Personal_Care Automotive Retail Professional_Services Entertainment Travel_Hospitality Education Home_Services Healthcare Veterinary_Services Sum_of_Industries
tUFrWirKiKi_TAnsVWINQQ Target 5255 E Broadway Blvd Tucson AZ 85711 32.22324 -110.8805 3.5 22 0 False True TRUE 2 FALSE False False FALSE u’no’ {‘garage’: False ‘street’: False ‘validated’: False ‘lot’: True ‘valet’: False} TRUE FALSE False FALSE False False NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA Department Stores, Shopping, Fashion, Home & Garden, Electronics, Furniture Stores 8:0-22:0 8:0-22:0 8:0-22:0 8:0-22:0 8:0-23:0 8:0-23:0 8:0-22:0 0 0 0 0 1 0 0 0 0,0,0,0,1
wghnIlMb_i5U46HMBGx9ig China Dragon Restaurant 1625 W Valencia Rd, Ste 101-103 Tucson AZ 85746 32.13230 -111.0000 3.0 23 0 NA True NA 2 NA True False FALSE ‘no’ {‘garage’: False ‘street’: False ‘validated’: False ‘lot’: False ‘valet’: False} NA NA False FALSE True NA ‘beer_and_wine’ True u’casual’ {‘romantic’: False ‘intimate’: False ‘classy’: False ‘hipster’: False ‘divey’: False ‘touristy’: False ‘trendy’: False ‘upscale’: False ‘casual’: True} NA True NA u’quiet’ NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA Restaurants, Chinese 11:0-21:0 11:0-21:0 11:0-21:0 11:0-21:0 11:0-21:0 11:0-21:0 11:0-21:0 1,0,0,0,0,0,0,0,0,0,0,0,1
txyXRytGjwOXvS8s4sc-WA Smoothie King 1070 E Tucson Marketplace Blvd Tucson AZ 85713 32.18679 -110.9548 3.0 29 1 NA True TRUE 2 NA True True NA u’free’ {u’valet’: False u’garage’: False u’street’: None u’lot’: True u’validated’: False} FALSE NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA Vitamins & Supplements, Ice Cream & Frozen Yogurt, Food, Juice Bars & Smoothies, Shopping 0:0-0:0 7:0-21:0 7:0-19:0 7:0-19:0 7:0-19:0 9:0-20:0 11:0-18:0 1 0 0 0 1 0 0 0 0,0,0,0,2
eYxGFkxo6m3SYGVTh5m2nQ Big Boyz Toyz Motorcycle Rentals 4158 E Grant Rd Tucson AZ 85712 32.25032 -110.9037 4.5 8 1 NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA Towing, Hotels & Travel, Automotive, Motorcycle Rental, Motorcycle Gear, Shopping 8:30-18:0 8:30-18:0 8:30-18:0 8:30-18:0 8:30-18:0 9:0-16:0 NA 0 0 0 1 1 0 0 1 0 0 0 0 3
x9K0RfZaT_zlw6DklBDzjw Gavi Italian Restaurant 7401 N La Cholla Blvd, Ste 146 Tucson AZ 85707 32.22167 -110.9258 3.5 9 0 NA True NA 2 NA True False NA NA {‘garage’: False ‘street’: False ‘validated’: False ‘lot’: True ‘valet’: False} NA NA True NA False NA u’full_bar’ True u’casual’ NA NA True NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA Italian, Restaurants NA NA NA NA NA NA NA 1 0 0 0 0 0 0 0 0,0,0,0,1
luGY5CJlb4CJ3LrBcSXMPA Desert Design Center 4146 E Speedway Blvd Tucson AZ 85712 32.23587 -110.9037 3.0 5 1 NA True NA 3 NA NA NA NA NA {‘garage’: False ‘street’: False ‘validated’: False ‘lot’: False ‘valet’: False} NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA Furniture Stores, Shopping, Home & Garden 9:30-18:0 9:30-18:0 9:30-18:0 9:30-18:0 9:30-18:0 9:30-18:0 12:0-17:0 0 0 0 0 1 0 0 0 0,0,0,0,1
wEfK2HIYNzkA1QYfUjHBaw Leslie’s 7955 E Broadway Blvd, Ste 147 Tucson AZ 85710 32.22154 -110.8221 2.5 15 1 NA True TRUE 2 NA NA NA NA NA None NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA Hardware Stores, Home Services, Hot Tub & Pool, Pool & Hot Tub Service, Home & Garden, Shopping 9:0-17:0 7:0-18:0 7:0-18:0 7:0-18:0 7:0-18:0 9:0-18:0 10:0-17:0 0 0 0 0 1 0 0 0 0 0 0 0 1
anLQj9AM8vjbcLSIE0iUgg Papa Murphy’s 7250 North La Cholla, Suite 186 Tucson AZ 85741 32.33857 -111.0108 4.0 22 1 NA True FALSE 1 NA True True NA NA {‘garage’: False ‘street’: False ‘validated’: False ‘lot’: False ‘valet’: False} NA NA False TRUE False NA u’none’ True u’casual’ NA NA True True u’quiet’ NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA Pizza, Restaurants NA NA NA NA NA NA NA 1 0 0 0 0 0 0 0 0,0,0,0,1
9C2rpb56aQvW0ViZHK9sPw Home Plate Sports Pub 4880 E 22nd St Tucson AZ 85711 32.20631 -110.8878 3.5 56 1 NA True TRUE 1 NA NA False FALSE u’no’ {u’valet’: False u’garage’: False u’street’: None u’lot’: None u’validated’: None} NA TRUE True TRUE False NA u’full_bar’ True u’casual’ {u’divey’: None u’hipster’: None u’casual’: True u’touristy’: None u’trendy’: None u’intimate’: False u’romantic’: None u’classy’: False u’upscale’: None} NA True NA u’average’ NA NA NA {u’dj’: None u’live’: False u’jukebox’: None u’video’: False u’background_music’: False u’karaoke’: None u’no_music’: False} NA NA {u’monday’: True u’tuesday’: False u’wednesday’: False u’thursday’: False u’friday’: False u’saturday’: False u’sunday’: False} NA NA NA NA NA,NA,NA,NA,“Bars, Restaurants, Sports Bars, Nightlife”,10:0-2:0,10:0-2:0,10:0-2:0,10:0-2:0,10:0-2:0,10:0-2:0,10:0-2:0,1,0,0,0,0,0,0,0,0,0,0,0,1
E_kLUlU60fzoZYAOxbH_3w Pain Institute of Southern Arizona 4582 N 1st Ave, Ste 170 Tucson AZ 85718 32.28918 -110.9572 2.5 43 1 True True NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA True NA NA NA NA NA NA NA NA NA Pain Management, Doctors, Health & Medical, Osteopathic Physicians 7:0-16:30 7:0-16:30 7:0-16:30 7:0-16:30 7:0-16:30 NA NA 0 0 0 0 0 0 0 0 0 0 1 0 1 
# Check review data structure  
cat("\nReview data sample:\n")
## 
## Review data sample:
knitr::kable(head(df_review, 10), 
             caption = "Sample of Review Data")
Sample of Review Data
review_id user_id business_id stars useful funny cool date
saUsX_uimxRlCVr67Z4Jig 8g_iMtfSiwikVnbP2etR0A YjUWPpI6HXG530lwP-fb2A 3 0 0 0 2014-02-05 20:30:30
UBp0zWyH60Hmw6Fsasei7w 4Uh27DgGzsp6PqrH913giQ otQS34_MymijPTdNBoBdCw 4 0 2 0 2011-10-27 17:12:05
ymhbOMW63B_vGaRFR3XT0A yZdAhQ_KmKuCZmbBy9YDQg 5Ce3lZksYVkCbrihqylVHQ 5 0 0 0 2014-07-25 17:56:26
onlgwy5qGDEzddsrnIvtWg pYXeL0RCqus2IfhthYCOyA W7NxQw8UYFR0HLPrI08tvw 4 0 0 0 2012-02-01 14:21:25
ypFqmURIY41F4pWaMW1VrQ hToW2eoTAYlOcuNz7I73Cg 7L1kXfwU5XM6f6Jxo_L7Yg 5 1 0 0 2015-07-03 22:36:04
_yfUDZFGYfc4Lmq2l-U7ZA jVVNbLCq2PaxMYtdNkxWGA VO5Jphyqmoa62lbA7-vzpQ 5 0 0 0 2013-05-07 15:19:23
M7h8K7PrLrN8yQb3RhHLdA CysTz612pTix75Fye6-D4Q WLiqfxv_GhFFA5sm878a2w 5 0 1 1 2013-12-29 03:15:34
dEVDg2FyOR2sJnbR_JcivQ 7zJVt_Bl1h5hLjhxl2c1PQ 3YERGr7UbpSpddqL0Eiu5g 5 1 0 1 2017-12-16 00:18:10
_u331rYIp7qJ7aWvLup9Ug BUpAhzsbbKbTqyBjT8YnNQ wa_bwyY57etHjtJ2Fw0E3g 5 0 0 0 2014-04-01 05:00:09
RgonT2qGpUte6RFpoSSxuQ iaLZ0oJ3sqCpnW-PY5_i4w dygflr6j51vf3-lM5pUIkw 3 0 0 0 2017-04-07 21:52:05 
# Check tip data structure  
cat("\nTip data sample:\n")
## 
## Tip data sample:
knitr::kable(head(df_tip, 10), 
             caption = "Sample of Tip Data")
Sample of Tip Data
user_id business_id date compliment_count
AGNUgVwnZUey3gcPCJ76iw 3uLgwr0qeCNMjKenHJwPGQ 2012-05-18 02:17:21 0
1viszljzsa8W7VnYPuWosw Zg83olSjsdXAhE5EEy5IcQ 2018-04-01 01:40:12 0
evVjF1TOaNwuUK7oY2k53w 5yjHCLzwHgjKNNHgAaUC8Q 2012-03-28 23:27:41 0
ztPinMY7A5o4HteSMyBzdw 1lSj24GO6VMJ2uB6hW3pGA 2017-03-17 21:35:00 0
juigu8XxwzLAeikeBzULbg DlOZoZZRLEOuh9FdAw1few 2013-08-18 23:54:19 0
8ohYIH9wHjKdWc4yAROqBw 3w7NRntdQ9h0KwDsksIt5Q 2013-03-09 19:23:16 0
td7JrqJMybMZKKoa1h9Qsg tV46IhCfHbsx_af-pMupiw 2012-05-30 20:19:18 0
wzX8UssOglKwSDx8qob8zA yZwlpIhzVRxVdbIkuw9-VA 2015-09-10 05:31:03 0
4Dy4PVb5lykv4N0kYiRhVA HpmsY3JJ4gIuS8EiVHwMJw 2017-11-18 21:11:26 0
3vofzamIMkrS0PPf5aSmRQ H5nMeiXJmLpBALeV8A-FeQ 2011-08-04 04:24:32 0

Feature Engineering from Review Data

# Create aggregated review features by business_id

review_features <- df_review |> 
  group_by(business_id) |>
  summarise(
    avg_review_stars = mean(stars, na.rm = TRUE),
    total_reviews = n(),
    avg_useful = mean(useful, na.rm = TRUE),
    avg_funny = mean(funny, na.rm = TRUE),
    avg_cool = mean(cool, na.rm = TRUE),
    most_recent_review = max(date, na.rm = TRUE), 
    .groups = 'drop'
  )

cat("Review features created for", nrow(review_features), "businesses\n")
## Review features created for 6187 businesses
knitr::kable(head(review_features, 10), 
             caption = "Aggregated Review Features by Business", 
             digits = 3)
Aggregated Review Features by Business
business_id avg_review_stars total_reviews avg_useful avg_funny avg_cool most_recent_review
–9osgUCSDUWUkoTLdvYhQ 4.800 30 0.800 0.133 0.633 2022-01-07 16:50:23
–MbOh2O1pATkXa7xbU6LA 3.963 27 0.778 0.259 0.259 2017-11-12 01:40:47
–Sd93OFWITqDHifMoAy9g 4.333 6 0.833 0.000 0.000 2019-06-05 22:44:05
-02xFuruu85XmDn2xiynJw 4.686 121 0.744 0.066 0.124 2021-06-24 15:59:50
-0E708CsutslDLUuuzWGYQ 2.208 24 2.000 0.333 0.167 2020-01-16 02:34:38
-0G_6-KFGpCpxTUlVXCMYQ 3.578 154 2.429 0.474 1.234 2015-03-07 23:46:33
-0yr__Vo4hOCZwvzRWd1zA 4.111 9 1.444 0.111 0.333 2019-07-13 21:11:15
-1MhPXk1FglglUAmuPLIGg 3.877 114 0.825 0.272 0.544 2021-11-30 15:15:15
-1iLbEf1NwY-OJp5Hg-3Sg 3.154 13 0.923 0.154 0.308 2019-10-23 01:35:34
-1w9JMktu9oWTXwNqtZQoA 3.882 17 0.824 0.118 0.412 2021-07-06 14:02:22

Merge Business and Review Data

# Merge business data with review features
df_merged <- df_business |> 
  left_join(review_features, by = "business_id") |>
  filter(!is.na(avg_review_stars)) # remove businesses without reviews

cat("Merged dataset shape:", dim(df_merged)[1], "x", dim(df_merged)[2], "\n")
## Merged dataset shape: 6187 x 77
# Check what columns we actually have
cat("\nColumns in merged dataset:\n")
## 
## Columns in merged dataset:
print(colnames(df_merged))
##  [1] "business_id"                          
##  [2] "name"                                 
##  [3] "address"                              
##  [4] "city"                                 
##  [5] "state"                                
##  [6] "postal_code"                          
##  [7] "latitude"                             
##  [8] "longitude"                            
##  [9] "stars"                                
## [10] "review_count"                         
## [11] "is_open"                              
## [12] "attributes.ByAppointmentOnly"         
## [13] "attributes.BusinessAcceptsCreditCards"
## [14] "attributes.BikeParking"               
## [15] "attributes.RestaurantsPriceRange2"    
## [16] "attributes.CoatCheck"                 
## [17] "attributes.RestaurantsTakeOut"        
## [18] "attributes.RestaurantsDelivery"       
## [19] "attributes.Caters"                    
## [20] "attributes.WiFi"                      
## [21] "attributes.BusinessParking"           
## [22] "attributes.WheelchairAccessible"      
## [23] "attributes.HappyHour"                 
## [24] "attributes.OutdoorSeating"            
## [25] "attributes.HasTV"                     
## [26] "attributes.RestaurantsReservations"   
## [27] "attributes.DogsAllowed"               
## [28] "attributes.Alcohol"                   
## [29] "attributes.GoodForKids"               
## [30] "attributes.RestaurantsAttire"         
## [31] "attributes.Ambience"                  
## [32] "attributes.RestaurantsTableService"   
## [33] "attributes.RestaurantsGoodForGroups"  
## [34] "attributes.DriveThru"                 
## [35] "attributes.NoiseLevel"                
## [36] "attributes.GoodForMeal"               
## [37] "attributes.BusinessAcceptsBitcoin"    
## [38] "attributes.Smoking"                   
## [39] "attributes.Music"                     
## [40] "attributes.GoodForDancing"            
## [41] "attributes.AcceptsInsurance"          
## [42] "attributes.BestNights"                
## [43] "attributes.BYOB"                      
## [44] "attributes.Corkage"                   
## [45] "attributes.BYOBCorkage"               
## [46] "attributes.HairSpecializesIn"         
## [47] "attributes.Open24Hours"               
## [48] "attributes.RestaurantsCounterService" 
## [49] "attributes.AgesAllowed"               
## [50] "attributes.DietaryRestrictions"       
## [51] "categories"                           
## [52] "hours.Monday"                         
## [53] "hours.Tuesday"                        
## [54] "hours.Wednesday"                      
## [55] "hours.Thursday"                       
## [56] "hours.Friday"                         
## [57] "hours.Saturday"                       
## [58] "hours.Sunday"                         
## [59] "Food_Beverage"                        
## [60] "Health_Fitness"                       
## [61] "Personal_Care"                        
## [62] "Automotive"                           
## [63] "Retail"                               
## [64] "Professional_Services"                
## [65] "Entertainment"                        
## [66] "Travel_Hospitality"                   
## [67] "Education"                            
## [68] "Home_Services"                        
## [69] "Healthcare"                           
## [70] "Veterinary_Services"                  
## [71] "Sum_of_Industries"                    
## [72] "avg_review_stars"                     
## [73] "total_reviews"                        
## [74] "avg_useful"                           
## [75] "avg_funny"                            
## [76] "avg_cool"                             
## [77] "most_recent_review"
# Check for missing values in key columns that exist
key_cols <- c("stars", "avg_review_stars", "total_reviews", "avg_useful", "avg_funny", "avg_cool")
existing_cols <- key_cols[key_cols %in% colnames(df_merged)]

if(length(existing_cols) > 0) {
  cat("\nMissing values in key columns:\n")
  
missing_summary <- df_merged[existing_cols] |>
  summarise_all(~sum(is.na(.)))
  
  knitr::kable(missing_summary, 
               caption = "Missing Values Summary", 
               digits = 0)
} else {
  cat("Key columns not found in dataset. Please check column names.\n")
}
## 
## Missing values in key columns:
Missing Values Summary
stars avg_review_stars total_reviews avg_useful avg_funny avg_cool
0 0 0 0 0 0

Target Variable Creation

# Create binary classification target (High vs. Low) rating
df_merged <- df_merged |>
  mutate(rating_category = case_when(    
    avg_review_stars < 3.0 ~ "Low",
    avg_review_stars >= 3.0 & stars < 4.0 ~ "Medium",
    avg_review_stars >= 4.0 ~ "High"  
),  
# Remove Medium category for binary classification   
binary_target = case_when (
  rating_category == "Low" ~ 0,
  rating_category == "High" ~ 1
  )
) |>
  filter(rating_category %in% c("Low", "High")) # Keep only Low and High

cat("Target variable distribution:\n")
## Target variable distribution:
target_dist <- table(df_merged$rating_category)
knitr::kable(as.data.frame(target_dist), 
             caption = "Rating Category Distribution",
             col.names = c("Rating Category", "Count"))
Rating Category Distribution
Rating Category Count
High 2337
Low 1547 
cat("\nBinary target distribution:\n")
## 
## Binary target distribution:
binary_dist <- table(df_merged$binary_target)
knitr::kable(as.data.frame(binary_dist), 
             caption = "Binary Target Distribution",
             col.names = c("Binary Target (0=Low, 1=High)", "Count"))
Binary Target Distribution
Binary Target (0=Low, 1=High) Count
0 1547
1 2337 
# Regression target is the numeric 'stars' column
cat("\nRegression target (stars) summary:\n")
## 
## Regression target (stars) summary:
stars_summary <- df_merged |> 
  summarise(
    Min = min(stars, na.rm = TRUE),
    Q1 = quantile(stars, 0.25, na.rm = TRUE),
    Median = median(stars, na.rm = TRUE),
    Mean = mean(stars, na.rm = TRUE),
    Q3 = quantile(stars, 0.75, na.rm = TRUE),
    Max = max(stars, na.rm = TRUE),
    SD = sd(stars, na.rm = TRUE)
  )

knitr::kable(stars_summary, 
             caption = "Summary Statistics for Business Star Ratings",
             digits = 3)
Summary Statistics for Business Star Ratings
Min Q1 Median Mean Q3 Max SD
1 2.5 4 3.63 4.5 5 1.117

Feature Selection and Preparation

# First, let's check what columns are actually in our merged dataset
cat("Columns in df_merged:\n")
## Columns in df_merged:
print(colnames(df_merged))
##  [1] "business_id"                          
##  [2] "name"                                 
##  [3] "address"                              
##  [4] "city"                                 
##  [5] "state"                                
##  [6] "postal_code"                          
##  [7] "latitude"                             
##  [8] "longitude"                            
##  [9] "stars"                                
## [10] "review_count"                         
## [11] "is_open"                              
## [12] "attributes.ByAppointmentOnly"         
## [13] "attributes.BusinessAcceptsCreditCards"
## [14] "attributes.BikeParking"               
## [15] "attributes.RestaurantsPriceRange2"    
## [16] "attributes.CoatCheck"                 
## [17] "attributes.RestaurantsTakeOut"        
## [18] "attributes.RestaurantsDelivery"       
## [19] "attributes.Caters"                    
## [20] "attributes.WiFi"                      
## [21] "attributes.BusinessParking"           
## [22] "attributes.WheelchairAccessible"      
## [23] "attributes.HappyHour"                 
## [24] "attributes.OutdoorSeating"            
## [25] "attributes.HasTV"                     
## [26] "attributes.RestaurantsReservations"   
## [27] "attributes.DogsAllowed"               
## [28] "attributes.Alcohol"                   
## [29] "attributes.GoodForKids"               
## [30] "attributes.RestaurantsAttire"         
## [31] "attributes.Ambience"                  
## [32] "attributes.RestaurantsTableService"   
## [33] "attributes.RestaurantsGoodForGroups"  
## [34] "attributes.DriveThru"                 
## [35] "attributes.NoiseLevel"                
## [36] "attributes.GoodForMeal"               
## [37] "attributes.BusinessAcceptsBitcoin"    
## [38] "attributes.Smoking"                   
## [39] "attributes.Music"                     
## [40] "attributes.GoodForDancing"            
## [41] "attributes.AcceptsInsurance"          
## [42] "attributes.BestNights"                
## [43] "attributes.BYOB"                      
## [44] "attributes.Corkage"                   
## [45] "attributes.BYOBCorkage"               
## [46] "attributes.HairSpecializesIn"         
## [47] "attributes.Open24Hours"               
## [48] "attributes.RestaurantsCounterService" 
## [49] "attributes.AgesAllowed"               
## [50] "attributes.DietaryRestrictions"       
## [51] "categories"                           
## [52] "hours.Monday"                         
## [53] "hours.Tuesday"                        
## [54] "hours.Wednesday"                      
## [55] "hours.Thursday"                       
## [56] "hours.Friday"                         
## [57] "hours.Saturday"                       
## [58] "hours.Sunday"                         
## [59] "Food_Beverage"                        
## [60] "Health_Fitness"                       
## [61] "Personal_Care"                        
## [62] "Automotive"                           
## [63] "Retail"                               
## [64] "Professional_Services"                
## [65] "Entertainment"                        
## [66] "Travel_Hospitality"                   
## [67] "Education"                            
## [68] "Home_Services"                        
## [69] "Healthcare"                           
## [70] "Veterinary_Services"                  
## [71] "Sum_of_Industries"                    
## [72] "avg_review_stars"                     
## [73] "total_reviews"                        
## [74] "avg_useful"                           
## [75] "avg_funny"                            
## [76] "avg_cool"                             
## [77] "most_recent_review"                   
## [78] "rating_category"                      
## [79] "binary_target"
# Check the structure of key industry columns to see what we're dealing with
cat("\nChecking industry column types:\n")
## 
## Checking industry column types:
industry_cols <- c("Food_Beverage", "Retail", "Healthcare", "Entertainment", "Health_Fitness")
existing_industry_cols <- industry_cols[industry_cols %in% colnames(df_merged)]

if(length(existing_industry_cols) > 0) {
  for(col in existing_industry_cols) {
    cat(paste("Column:", col, "- Class:", class(df_merged[[col]]), "- Sample values:", 
              paste(head(unique(df_merged[[col]]), 3), collapse = ", "), "\n"))
  }
}
## Column: Food_Beverage - Class: character - Sample values: 0, 9:30-18:0, 9:0-17:0 
## Column: Retail - Class: character - Sample values: 1, 0, NA 
## Column: Healthcare - Class: character - Sample values: 0, 1, NA 
## Column: Entertainment - Class: character - Sample values: 0, NA, 1 
## Column: Health_Fitness - Class: character - Sample values: 0, 9:30-18:0, 10:0-15:0
# Select predictor variables that we know exist and are properly formatted
# Start with basic business characteristics that should be numeric
safe_predictor_cols <- c()

# Add review count if it exists
if("review_count" %in% colnames(df_merged)) {
  safe_predictor_cols <- c(safe_predictor_cols, "review_count")
}

# Add location data if it exists  
if("latitude" %in% colnames(df_merged) && "longitude" %in% colnames(df_merged)) {
  safe_predictor_cols <- c(safe_predictor_cols, "latitude", "longitude")
}

# Add aggregated review features (but exclude avg_review_stars to avoid data leakage)
review_feature_cols <- c("total_reviews", "avg_useful", "avg_funny", "avg_cool")
existing_review_cols <- review_feature_cols[review_feature_cols %in% colnames(df_merged)]
safe_predictor_cols <- c(safe_predictor_cols, existing_review_cols)

# Only add industry columns if they're properly formatted (numeric 0/1)
for(col in existing_industry_cols) {
  if(is.numeric(df_merged[[col]]) && all(df_merged[[col]] %in% c(0, 1, NA))) {
    safe_predictor_cols <- c(safe_predictor_cols, col)
    cat("Added industry column:", col, "\n")
  } else {
    cat("Skipping problematic industry column:", col, "\n")
  }
}
## Skipping problematic industry column: Food_Beverage 
## Skipping problematic industry column: Retail 
## Skipping problematic industry column: Healthcare 
## Skipping problematic industry column: Entertainment 
## Skipping problematic industry column: Health_Fitness
cat("\nUsing predictor columns:", paste(safe_predictor_cols, collapse = ", "), "\n")
## 
## Using predictor columns: review_count, latitude, longitude, total_reviews, avg_useful, avg_funny, avg_cool
# Create feature matrix with only safe, properly formatted columns
X <- df_merged[safe_predictor_cols] |>
  mutate_all(~as.numeric(as.character(.))) |>  # Ensure all are numeric  
  mutate_all(~ifelse(is.na(.), 0, .))  # Replace NAs with 0# Target variables
y_binary <- df_merged$binary_target
y_regression <- df_merged$avg_review_stars

cat("Feature matrix shape:", dim(X)[1], "x", dim(X)[2], "\n")
## Feature matrix shape: 3884 x 7
cat("Features used:\n")
## Features used:
# Create a nice table showing the features we're actually using
if(ncol(X) > 0) {
  feature_info <- data.frame(
    Feature_Name = colnames(X),
    Feature_Type = ifelse(colnames(X) %in% existing_industry_cols, "Industry Indicator",
                         ifelse(colnames(X) %in% c("review_count", "latitude", "longitude"), "Business Characteristic",
                               "Aggregated Review Metric")),
    Data_Type = sapply(X, class),
    Sample_Range = sapply(X, function(x) paste(round(range(x, na.rm = TRUE), 2), collapse = " to "))
  )
  
  knitr::kable(feature_info, 
               caption = "Predictor Variables Used in Models",
               col.names = c("Feature Name", "Type", "Data Type", "Range"))
} else {
  cat("No valid predictor columns found!\n")
}
Predictor Variables Used in Models
Feature Name Type Data Type Range
review_count review_count Business Characteristic numeric 5 to 2126
latitude latitude Business Characteristic numeric 31.9 to 32.53
longitude longitude Business Characteristic numeric -111.22 to -110.63
total_reviews total_reviews Aggregated Review Metric numeric 5 to 2230
avg_useful avg_useful Aggregated Review Metric numeric 0 to 16.55
avg_funny avg_funny Aggregated Review Metric numeric 0 to 5.27
avg_cool avg_cool Aggregated Review Metric numeric 0 to 5.82

Data Splitting

# Split data for classification
train_idx_cls <- createDataPartition(y_binary, p = 0.8, list = FALSE)
X_train_cls <- X[train_idx_cls, ]
X_test_cls <- X[-train_idx_cls, ]
y_train_cls <- y_binary[train_idx_cls]
y_test_cls <- y_binary[-train_idx_cls]

# Split data for regression  
train_idx_reg <- createDataPartition(y_regression, p = 0.8, list = FALSE)
X_train_reg <- X[train_idx_reg, ]
X_test_reg <- X[-train_idx_reg, ]
y_train_reg <- y_regression[train_idx_reg]
y_test_reg <- y_regression[-train_idx_reg]

cat("Classification - Training set:", dim(X_train_cls)[1], "Test set:", dim(X_test_cls)[1], "\n")
## Classification - Training set: 3108 Test set: 776
cat("Regression - Training set:", dim(X_train_reg)[1], "Test set:", dim(X_test_reg)[1], "\n")
## Regression - Training set: 3110 Test set: 774

Binary Classification Model

cat("=" %>% rep(60) %>% paste(collapse = ""), "\n")
## ============================================================
cat("BINARY CLASSIFICATION: High vs Low Rating\n")
## BINARY CLASSIFICATION: High vs Low Rating
cat("=" %>% rep(60) %>% paste(collapse = ""), "\n")
## ============================================================
# Prepare training data - convert binary target to proper factor levels
train_data_cls <- cbind(X_train_cls, target = factor(y_train_cls, levels = c(0, 1), labels = c("Low", "High")))

# Set up cross-validation
ctrl <- trainControl(method = "cv", number = 5, classProbs = TRUE, summaryFunction = twoClassSummary)

# Train Random Forest model
rf_model_cls <- train(target ~ ., 
                      data = train_data_cls,
                      method = "rf",
                      trControl = ctrl,
                      tuneGrid = expand.grid(mtry = c(2, 3, 4, 5)),
                      ntree = 100,
                      metric = "ROC")

cat("Best parameters:\n")
## Best parameters:
print(rf_model_cls$bestTune)
##   mtry
## 3    4
cat("Best cross-validation ROC:", round(max(rf_model_cls$results$ROC), 4), "\n")
## Best cross-validation ROC: 0.8604
# Predict on test set
y_pred_cls <- predict(rf_model_cls, newdata = X_test_cls)
y_pred_prob <- predict(rf_model_cls, newdata = X_test_cls, type = "prob")[,"High"]

# Convert test labels to same factor format for evaluation
y_test_cls_factor <- factor(y_test_cls, levels = c(0, 1), labels = c("Low", "High"))

# Evaluate model
conf_matrix <- confusionMatrix(y_pred_cls, y_test_cls_factor, positive = "High")
print(conf_matrix)
## Confusion Matrix and Statistics
## 
##           Reference
## Prediction Low High
##       Low  217   90
##       High  83  386
##                                              
##                Accuracy : 0.7771             
##                  95% CI : (0.7461, 0.8059)   
##     No Information Rate : 0.6134             
##     P-Value [Acc > NIR] : <0.0000000000000002
##                                              
##                   Kappa : 0.532              
##                                              
##  Mcnemar's Test P-Value : 0.6483             
##                                              
##             Sensitivity : 0.8109             
##             Specificity : 0.7233             
##          Pos Pred Value : 0.8230             
##          Neg Pred Value : 0.7068             
##              Prevalence : 0.6134             
##          Detection Rate : 0.4974             
##    Detection Prevalence : 0.6044             
##       Balanced Accuracy : 0.7671             
##                                              
##        'Positive' Class : High               
## 
# Calculate additional metrics
precision <- conf_matrix$byClass['Precision']
recall <- conf_matrix$byClass['Recall'] 
f1_score <- conf_matrix$byClass['F1']
accuracy <- conf_matrix$overall['Accuracy']

cat("\nClassification Metrics:\n")
## 
## Classification Metrics:
cat("Accuracy:", round(accuracy, 4), "\n")
## Accuracy: 0.7771
cat("Precision:", round(precision, 4), "\n")
## Precision: 0.823
cat("Recall:", round(recall, 4), "\n")
## Recall: 0.8109
cat("F1 Score:", round(f1_score, 4), "\n")
## F1 Score: 0.8169

Regression Model

cat("=" %>% rep(60) %>% paste(collapse = ""), "\n")
## ============================================================
cat("REGRESSION: Predicting Star Ratings\n") 
## REGRESSION: Predicting Star Ratings
cat("=" %>% rep(60) %>% paste(collapse = ""), "\n")
## ============================================================
# Prepare training data
train_data_reg <- cbind(X_train_reg, target = y_train_reg)

# Set up cross-validation for regression
ctrl_reg <- trainControl(method = "cv", number = 5)

# Train Gradient Boosting model
gbm_model_reg <- train(target ~ ., 
                       data = train_data_reg,
                       method = "gbm",
                       trControl = ctrl_reg,
                       tuneGrid = expand.grid(
                         interaction.depth = c(3, 5),
                         n.trees = c(100, 200),
                         shrinkage = c(0.01, 0.1),
                         n.minobsinnode = 10
                       ),
                       metric = "RMSE",
                       verbose = FALSE)

cat("Best parameters:\n")
## Best parameters:
print(gbm_model_reg$bestTune)
##   n.trees interaction.depth shrinkage n.minobsinnode
## 5     100                 3       0.1             10
cat("Best cross-validation RMSE:", round(min(gbm_model_reg$results$RMSE), 4), "\n")
## Best cross-validation RMSE: 0.8668
# Predict on test set
y_pred_reg <- predict(gbm_model_reg, newdata = X_test_reg)

# Evaluate regression model
mse <- mean((y_test_reg - y_pred_reg)^2)
rmse <- sqrt(mse)
mae <- mean(abs(y_test_reg - y_pred_reg))
r2 <- cor(y_test_reg, y_pred_reg)^2

cat("\nRegression Metrics:\n")
## 
## Regression Metrics:
cat("MSE:", round(mse, 4), "\n")
## MSE: 0.7227
cat("RMSE:", round(rmse, 4), "\n")
## RMSE: 0.8501
cat("MAE:", round(mae, 4), "\n")
## MAE: 0.6675
cat("R² Score:", round(r2, 4), "\n")
## R² Score: 0.3983

Feature Importance Analysis

# Classification feature importance
importance_cls <- varImp(rf_model_cls)$importance
importance_df_cls <- data.frame(
  feature = rownames(importance_cls),
  importance = importance_cls$Overall
) |>
  arrange(desc(importance))

cat("Top 5 Important Features for Classification:\n")
## Top 5 Important Features for Classification:
knitr::kable(head(importance_df_cls, 10), 
             caption = "Feature Importance for Classification Model",
             col.names = c("Feature", "Importance"),
             digits = 3, row.names = FALSE)
Feature Importance for Classification Model
Feature Importance
avg_cool 100.000
avg_funny 48.724
latitude 26.260
longitude 22.586
avg_useful 21.392
total_reviews 0.638
review_count 0.000 
# Regression feature importance  
importance_reg <- varImp(gbm_model_reg)$importance
importance_df_reg <- data.frame(
  feature = rownames(importance_reg),
  importance = importance_reg$Overall
) |>
  arrange(desc(importance))

cat("\nTop 5 Important Features for Regression:\n")
## 
## Top 5 Important Features for Regression:
knitr::kable(head(importance_df_reg, 10), 
             caption = "Feature Importance for Regression Model",
             col.names = c("Feature", "Importance"),
             digits = 3, row.names = FALSE)
Feature Importance for Regression Model
Feature Importance
avg_cool 100.000
avg_funny 70.323
avg_useful 14.590
latitude 11.157
longitude 7.810
review_count 0.310
total_reviews 0.000

Visualizations

# Create plots
par(mfrow = c(2, 3))

# 1. Confusion Matrix Heatmap
cm_matrix <- as.matrix(conf_matrix$table)
image(1:ncol(cm_matrix), 1:nrow(cm_matrix), t(cm_matrix[nrow(cm_matrix):1,]), 
      col = heat.colors(100), xlab = "Predicted", ylab = "Actual",
      main = "Confusion Matrix", axes = FALSE)
axis(1, at = 1:ncol(cm_matrix), labels = colnames(cm_matrix))
axis(2, at = 1:nrow(cm_matrix), labels = rev(rownames(cm_matrix)))
text(expand.grid(1:ncol(cm_matrix), 1:nrow(cm_matrix)), 
     sprintf("%d", as.vector(t(cm_matrix[nrow(cm_matrix):1,]))), cex = 1.2)

# 2. Feature Importance - Classification
barplot(head(importance_df_cls$importance, 8), 
        names.arg = head(importance_df_cls$feature, 8),
        main = "Top Features - Classification",
        ylab = "Importance", las = 2, cex.names = 0.8,
        col = "lightblue")

# 3. Feature Importance - Regression
barplot(head(importance_df_reg$importance, 8),
        names.arg = head(importance_df_reg$feature, 8), 
        main = "Top Features - Regression",
        ylab = "Importance", las = 2, cex.names = 0.8,
        col = "lightgreen")

# 4. Actual vs Predicted (Regression)
# Get the actual average review stars for the test set businesses
test_indices_reg <- train_idx_reg  # This should be the indices used for regression split
actual_avg_review_stars <- df_merged$avg_review_stars[-test_indices_reg]  # Test set actual review averages

plot(actual_avg_review_stars, y_pred_reg,
     main = paste("Predicted Business Stars vs Actual Review Average\nR² =", round(r2, 3)),
     xlab = "Actual Average Review Stars", ylab = "Predicted Business Stars",
     col = "darkblue", pch = 16, alpha = 0.6)
abline(a = 0, b = 1, col = "red", lwd = 2, lty = 2)
grid()

# 5. Distribution of Predictions
hist(y_pred_reg, breaks = 20, 
     main = "Distribution of Predicted Ratings",
     xlab = "Predicted Stars", col = "lightcoral", border = "black")

# 6. Residuals Plot
# Calculate residuals as: Actual Average Review Stars - Predicted Business Stars
residuals <- actual_avg_review_stars - y_pred_reg
plot(y_pred_reg, residuals,
     main = "Residuals vs Predicted Business Stars",
     xlab = "Predicted Business Stars", ylab = "Residuals (Actual Review Avg - Predicted)",
     col = "purple", pch = 16)
abline(h = 0, col = "red", lwd = 2, lty = 2)
grid()

par(mfrow = c(1, 1))

Results Summary

cat("=" %>% rep(60) %>% paste(collapse = ""), "\n")
## ============================================================
cat("FINAL RESULTS SUMMARY\n")
## FINAL RESULTS SUMMARY
cat("=" %>% rep(60) %>% paste(collapse = ""), "\n")
## ============================================================
cat("\nCLASSIFICATION RESULTS (High vs Low Rating):\n")
## 
## CLASSIFICATION RESULTS (High vs Low Rating):
cat("-" %>% rep(40) %>% paste(collapse = ""), "\n")
## ----------------------------------------
cat("Model: Random Forest\n")
## Model: Random Forest
cat("Accuracy:", round(accuracy, 4), "\n")
## Accuracy: 0.7771
cat("Precision:", round(precision, 4), "\n") 
## Precision: 0.823
cat("Recall:", round(recall, 4), "\n")
## Recall: 0.8109
cat("F1 Score:", round(f1_score, 4), "\n")
## F1 Score: 0.8169
cat("\nREGRESSION RESULTS (Star Rating Prediction):\n")
## 
## REGRESSION RESULTS (Star Rating Prediction):
cat("-" %>% rep(40) %>% paste(collapse = ""), "\n")
## ----------------------------------------
cat("Model: Gradient Boosting\n")
## Model: Gradient Boosting
cat("R² Score:", round(r2, 4), "\n")
## R² Score: 0.3983
cat("RMSE:", round(rmse, 4), "\n")
## RMSE: 0.8501
cat("MAE:", round(mae, 4), "\n")
## MAE: 0.6675
cat("\nKEY INSIGHTS:\n")
## 
## KEY INSIGHTS:
cat("-" %>% rep(15) %>% paste(collapse = ""), "\n")
## ---------------
cat("• Most important features include review-based metrics (avg_review_stars, total_reviews)\n")
## • Most important features include review-based metrics (avg_review_stars, total_reviews)
cat("• Industry type and location also play important roles in predictions\n")
## • Industry type and location also play important roles in predictions
cat("• Models show good performance for predicting business ratings from metadata\n")
## • Models show good performance for predicting business ratings from metadata
cat("• Review aggregation features significantly improve prediction accuracy\n")
## • Review aggregation features significantly improve prediction accuracy