# Suggestions to Limit Concussions on Punts For my analsis, I use logistic regression to estimate variables that most impact whether a punt is returned. I then propose small restrictions on the return team formations that will limit punt returns and estimate it will reduce concussions by approximately 29%. ## 1. Exploratory Analysis ```{r include=F, eval=T , echo=F} library(plyr) library(dplyr) library(data.table) library(ggplot2) library(ggthemes) library(gridExtra) library(stringr) library(fasttime) player_punt_data <- fread('../input/NFL-Punt-Analytics-Competition/player_punt_data.csv') play_information <- fread('../input/NFL-Punt-Analytics-Competition/play_information.csv') video_footage_injury <- fread('../input/NFL-Punt-Analytics-Competition/video_footage-injury.csv') game_data <- fread('../input/NFL-Punt-Analytics-Competition/game_data.csv') play_player_role_data <- fread('../input/NFL-Punt-Analytics-Competition/play_player_role_data.csv') player_punt_data <- fread('../input/NFL-Punt-Analytics-Competition/player_punt_data.csv') video_footage_control <- fread('../input/NFL-Punt-Analytics-Competition/video_footage-control.csv') video_review <- fread('../input/NFL-Punt-Analytics-Competition/video_review.csv') ``` The data contains 37 concussions from 6681 punt plays. The small sample size of concussions makes it hard to make any predictive models, but you can look at descriptive statistics/plots: ```{r include=F, eval=T , echo=F} game_data$concussion <- (game_data$GameKey %in% video_review$GameKey) game_data$concussion <- (game_data$GameKey %in% video_review$GameKey) game_data$Turf<-tolower(game_data$Turf) game_data$Turf[grepl("natural|natrual", game_data$Turf)]<-'grass' game_data$Turf[grepl("ubu ", game_data$Turf)]<-'ubu s5-m' game_data$Turf[grepl("fieldturf|field turf", game_data$Turf)]<-'fieldturf' game_data$Turf[grepl("artific|synthetic|turf", game_data$Turf)]<-'artificial' game_data[,.('count'=length(GameKey),'percent_concussion'=mean(concussion), 'count_concussion'=sum(concussion)), by=Turf][order(-percent_concussion)] game_data$Start_Time_Numeric<-as.numeric(sapply(strsplit(game_data$Start_Time, ":"), `[[`, 1))+ as.numeric(sapply(strsplit(game_data$Start_Time, ":"), `[[`, 2))/60 game_data$Start_Time_Bin<-cut(game_data$Start_Time_Numeric, breaks = c(10, 15,18, 24)) video_review$game_play_id <- paste(video_review$GameKey, video_review$PlayID, sep="_") play_information$game_play_id <- paste(play_information$GameKey, play_information$PlayID, sep="_") play_information$concussion <- (play_information$game_play_id %in% video_review$game_play_id) play_information[,.('count'=length(GameKey),'percent_concussion'=mean(concussion), 'count_concussion'=sum(concussion)), by=Quarter][order(-percent_concussion)] table(play_information$concussion) play_information$Game_Date<-as.Date(play_information$Game_Date, format="%m/%d/%") play_information$Punt.Success<-as.numeric(grepl("punts", play_information$PlayDescription)) play_information$Punt.Dist[play_information$Punt.Success==1]<-sapply(gregexpr("punts ", play_information$PlayDescription[play_information$Punt.Success==1]),`[[`, 1) play_information$Punt.Dist<-as.numeric(play_information$Punt.Dist) play_information$Punt.Dist<-as.numeric(substring(play_information$PlayDescription, play_information$Punt.Dist+6, play_information$Punt.Dist+7)) play_information$Returned<-as.numeric(grepl(" for ", play_information$PlayDescription)& grepl("punts" , play_information$PlayDescription)) play_information$Blocked<-as.numeric(grepl("BLOCKED" , play_information$PlayDescription)) play_information$Muffed<-as.numeric(grepl(" for ", play_information$PlayDescription)& grepl("punts" , play_information$PlayDescription)& grepl("MUFF", play_information$PlayDescription)) play_information$Return.Dist[play_information$Returned==1& play_information$Muffed==0]<-sapply(gregexpr("for ", play_information$PlayDescription[play_information$Returned==1& play_information$Muffed==0]),`[[`, 1) play_information$Return.Dist<-as.numeric(substring(play_information$PlayDescription, play_information$Return.Dist+4, play_information$Return.Dist+5)) play_information$Fair.Catch<-as.numeric(grepl("fair", play_information$PlayDescription)& grepl("catch" , play_information$PlayDescription)) play_information$Penalty<-play_information$PlayDescription play_information$Penalty[!grepl("PENALTY", play_information$Penalty) | grepl("Holding", play_information$Penalty)]<-NA play_information$Penalty[grepl("PENALTY", play_information$Penalty)]<- sapply(strsplit(play_information$Penalty[grepl("PENALTY", play_information$Penalty)], "PENALTY"), `[[`,2) play_information$Penalty[!is.na(play_information$Penalty)]<-sapply(strsplit(play_information$Penalty[!is.na(play_information$Penalty)], ", "), `[[`, 2) play_information$Dangerous.Penalty<-play_information$Penalty%in% c( "Illegal Block Above the Waist","Fair Catch Interference" , "Unnecessary Roughness", "Clipping", "Illegal Blindside Block" ) #reshape return type variable play_information$Return.Type<-ifelse(play_information$Muffed==1, "Muffed", ifelse(play_information$Returned==1, "Returned", ifelse(play_information$Fair.Catch==1, "Fair Catch", ifelse(play_information$Punt.Success==1, "Touchback/Out", "Other (Blocked/Penalty/Fumble)" )))) ###organize data ##### video_footage_injury$game_play_id <- paste(video_footage_injury$gamekey, video_footage_injury$playid, sep="_") video_review$Returned<-video_review$game_play_id%in% play_information$game_play_id[play_information$Returned==1] video_review$Primary_Number<-player_punt_data$Number[match(video_review$GSISID, player_punt_data$GSISID)] video_review$Partner_Number<-player_punt_data$Number[match(video_review$Primary_Partner_GSISID, player_punt_data$GSISID)] video_review$Penalty<-play_information$PlayDescription[match( video_review$game_play_id, play_information$game_play_id)] video_review$Penalty[!grepl("PENALTY", video_review$Penalty) | grepl("Holding", video_review$Penalty)]<-NA video_review$Penalty[grepl("PENALTY", video_review$Penalty)]<- sapply(strsplit(video_review$Penalty[grepl("PENALTY", video_review$Penalty)], "PENALTY"), `[[`,2) video_review$`PREVIEW LINK (5000K)`<-video_footage_injury$`PREVIEW LINK (5000K)`[match(video_footage_injury$game_play_id, video_review$game_play_id)] video_review[video_review$Primary_Impact_Type=="Helmet-to-helmet",] #player role data--has all player roles--can get formations by aggregating player roles play_player_role_data$game_play_id <- paste(play_player_role_data$GameKey, play_player_role_data$PlayID, sep="_") play_player_role_data <- play_player_role_data[order(game_play_id, Role)] play_formations <- play_player_role_data[,paste(Role, collapse="_"),by=game_play_id] colnames(play_formations) <- c('game_play_id','formation') head(play_formations) #full formation: too many combos to make much sense play_information <- merge(play_information, play_formations, by='game_play_id', all.x=T) play_player_role_data$punting_team <- (play_player_role_data$Role %in% c('GL','GR', "GLi","GLo", "GRi", "GRo",'PLW','PLT','PLG', 'PLS','PRG','PRT','PRW','PC','PPR','PPRi','PPRo','P', 'PPLi', 'PPLo','PPL')) play_formations_detail <- play_player_role_data[,paste(Role, collapse="_"),by=c('game_play_id', 'punting_team')] play_formations_detail<-reshape(play_formations_detail,timevar="punting_team",idvar="game_play_id",direction="wide") colnames(play_formations_detail) <- c('game_play_id','punt.team.formation', 'return.team.formation') play_information <- merge(play_information, play_formations_detail, by='game_play_id', all.x=T) video_review<-data.frame(video_review) video_review$concussion<-T all_player_roles<-merge(play_player_role_data, video_review[, c("GSISID", "concussion","game_play_id")], by=c("GSISID", "game_play_id"), all.x=T) all_player_roles<-data.frame(all_player_roles) all_player_roles$concussion[is.na(all_player_roles$concussion)]<-F video_review$caused.concussion<-T all_player_roles<-merge(all_player_roles, video_review[, c("Primary_Partner_GSISID", "caused.concussion","game_play_id")], by.x=c("GSISID", "game_play_id"),by.y=c("Primary_Partner_GSISID", "game_play_id"), all.x=T) all_player_roles$caused.concussion[is.na(all_player_roles$caused.concussion)]<-F all_player_roles$all.players<-T all_player_roles$all.players.video<-all_player_roles$game_play_id%in% video_review$game_play_id all_player_roles$Returned<-all_player_roles$game_play_id%in% play_information$game_play_id[play_information$Returned==1] all_player_roles<-data.frame(all_player_roles) all_player_roles$Position<-player_punt_data$Position[match(all_player_roles$GSISID, player_punt_data$GSISID)] all_player_roles$Role2<-gsub("1|2|3|4|5|6|7|8|9","", all_player_roles$Role) all_player_roles$Role3<-ifelse(grepl("PDR|PDL|PDM", all_player_roles$Role),"PD" , ifelse(grepl("PLR|PLM|PLL", all_player_roles$Role), "PL", ifelse(grepl("PC", all_player_roles$Role), "PC", ifelse(grepl("PLW|PRW", all_player_roles$Role), "PW", ifelse(grepl("PRG|PRT|PLG|PLT|PLS", all_player_roles$Role), "POL", ifelse(grepl("GL|GR",all_player_roles$Role), "G", ifelse(grepl("VR|VL", all_player_roles$Role), "V", all_player_roles$Role))))))) play_information$jammers<-str_count(play_information$return.team.formation, "V") play_information$PL<-str_count(play_information$return.team.formation, "PLL|PLM|PLR") play_information$PD<-str_count(play_information$return.team.formation, "PDL|PDR|PDM") play_information$FB<-str_count(play_information$return.team.formation, "PFB") play_information$Side.Overload<-abs(str_count(play_information$return.team.formation, "PDL")-str_count(play_information$return.team.formation, "PDR")+ str_count(play_information$return.team.formation, "PLL")-str_count(play_information$return.team.formation, "PLR") ) play_information$Box.Defenders<-play_information$PD+play_information$PL play_information$gunners<-str_count(play_information$punt.team.formation, "GL")+str_count(play_information$punt.team.formation, "GR") play_information$PP<-str_count(play_information$punt.team.formation, "PP") play_information$Box.Offense<-10-play_information$gunners-play_information$PP ``` ```{r include=T, eval=T , echo=F, fig.width=10, fig.height=4, fig.align="center"} ggplot(video_review, aes(x=Player_Activity_Derived, fill = Primary_Impact_Type))+ geom_bar(stat = 'count')+ ggtitle("Contact Type Resulting in Concussion")+ xlab(NULL) grouping.var<-"Role3" df.tidy<-melt(all_player_roles[, c("Role3", "Role2", "Role", "Position", "punting_team", "all.players", "all.players.video", "caused.concussion", "concussion")], id.vars=c("Role3", "Role2", "Role", "Position", "punting_team") ) df.tidy<-df.tidy[df.tidy$value==T,] df.tidy<-table(df.tidy$variable, df.tidy[, grouping.var]) counts<-melt(df.tidy) df.tidy<-melt(df.tidy/rowSums(df.tidy)) colnames(df.tidy)<-c("Player.Type", "Role", "Percent") # df.tidy<-df.tidy[df.tidy$Percent>.015,] #plot cdf of player type2 grouped by concussed, caused.concussion, and all.players ggplot(df.tidy, aes(Role, Percent, fill=Player.Type)) + geom_bar(position = "dodge", stat="identity") + labs(x = grouping.var, y="Percent of Player.Type with Each Role") + geom_text(aes(label=as.factor(counts$value)), vjust = -0.5, position = position_dodge(0.9)) ``` Helmet-to-helmet hits are now trying to be restricted, and they had been making up a good percentage of concussions. In addition, certain positions like players on the offensive line are getting concussed at a disproportionate amount. The small sample size, however, leads to lack of significance in much of the cross tabulations. For example, despite Preseason having a higher rate of concussions, a t-test shows how it's not statistically significant. ```{r include=T, eval=T , echo=T} t.test(play_information$concussion[which(play_information$Season_Type=="Pre")], play_information$concussion[which(play_information$Season_Type=="Reg")]) ``` One thing that is significant is that returned punts have a higher concussion rate, as is expected. ```{r include=T, eval=T , echo=T} t.test(play_information$concussion[which( play_information$Returned==1)], play_information$concussion[which( play_information$Returned==0)]) ``` ```{r include=F, eval=T , echo=F} ####NGS ANALYSIS---PLAYER IMPACT VARIABLE##### ngs.summarized<-data.frame(fread('../input/nfl-punt-data/NGS-summarized.csv')) ngs.summarized$max.impact[is.infinite(ngs.summarized$max.impact)]<-NA all_player_roles<-merge(data.frame(all_player_roles)[, !grepl("velocity|impact", colnames(all_player_roles))], ngs.summarized[, !colnames(ngs.summarized)%in%c("GameKey", "PlayID")], by=c("GSISID", "game_play_id"), all.x=T) all_player_roles<-data.table(all_player_roles) ngs.summarized.plays<-data.table(ngs.summarized) ngs.summarized.plays<-ngs.summarized.plays[, list(max.impact=max(max.impact, na.rm=T) ), by=c("game_play_id")] play_information<-merge(play_information, ngs.summarized.plays, by='game_play_id' ) ####NGS ANALYSIS-- DISTANCE TO LINE OF SCRIMMAGE##### #get line of scrimmage: ngs.organized<-fread('../input/nfl-punt-data/NGS-organized.csv') ngs.organized<-ngs.organized[order(ngs.organized$time_numeric, decreasing = F),] ngs.organized<-ngs.organized[ngs.organized$Event=='ball_snap',] ngs.organized<-ngs.organized[!duplicated(ngs.organized[, c("GSISID", "game_play_id")]),] cols<-!colnames(all_player_roles)%in% c("x", "y") all_player_roles<-merge(all_player_roles[,cols ,with=F], ngs.organized[, c("GSISID","x", "y", "game_play_id"), with=F], by=c("GSISID","game_play_id"), all.x=T) #ball.loc.x=location of long-snapper when ball is being snapped all_player_roles<-all_player_roles[, `:=`('ball.loc.y'=y[Role=="PLS"], 'ball.loc.x'=x[Role=="PLS"], 'p.loc.x'=x[Role=="P"] ), by=c("game_play_id")] #calculate distance to line of scrimmage ngs.organized<-fread('../input/nfl-punt-data/NGS-organized.csv') ngs.organized<-ngs.organized[order(ngs.organized$time_numeric, decreasing = F),] ngs.organized<-ngs.organized[ngs.organized$Event=='punt'& ngs.organized$TimeFromSnap>=0& ngs.organized$TimeFromSnap<=6,] ngs.organized<-ngs.organized[!duplicated(ngs.organized[, c("GSISID", "game_play_id")]),] all_player_roles<-merge(all_player_roles[, !colnames(all_player_roles)%in% c("x", "y"),with=F], ngs.organized[, c("GSISID","x", "y", "game_play_id"), with=F], by=c("GSISID","game_play_id"), all.x=T) all_player_roles$ball.loc.y[all_player_roles$ball.loc.y<=10]<-NA all_player_roles$dist.to.ball.y<-abs(all_player_roles$y-all_player_roles$ball.loc.y) all_player_roles$dist.to.ball.x<-ifelse( all_player_roles$p.loc.x=5), 'ball.loc.y'=abs(ball.loc.y[1]-26.65), 'Box.Defender.dist.at.punt.x.max'=max(dist.to.ball.x[!grepl("V|PFB", Role)& !Role=="PR"& punting_team==F], na.rm=T), 'PL.dist.at.punt.x'=mean(dist.to.ball.x[grepl("PLL|PLR|PLM", Role)& punting_team==F], na.rm=T), 'PL.dist.at.punt.x.med'=median(dist.to.ball.x[grepl("PLL|PLR|PLM", Role)& punting_team==F], na.rm=T), 'PL.dist.at.punt.x.max'=max(dist.to.ball.x[grepl("PLL|PLR|PLM", Role)& punting_team==F], na.rm=T) ), by="game_play_id"] players.downfield$Box.Defender.dist.at.punt.x.max[is.infinite(players.downfield$Box.Defender.dist.at.punt.x.max)]<-NA players.downfield$PL.dist.at.punt.x.max[is.infinite(players.downfield$PL.dist.at.punt.x.max)]<-NA play_information$Num.Downfield<-NA play_information<-merge(play_information[, !grepl("Downfield|dist[.]to|dist[.]from|at[.]punt|ball[.]loc", colnames(play_information)),with=F], players.downfield, by="game_play_id", all.x=T) ####NGS ANALYSIS-- DSITANCE TO PUNT RETURNER##### # must run "misc ngs organize.R" to get NGS csv file used here list_of_files<-list.files("../input/NFL-Punt-Analytics-Competition/",full.names = T) list_of_files<-list_of_files[grepl("NGS", list_of_files)] # conn<-unz("all.zip", file) # taconn[1] extractNGS<-function(file){ print(file) ngs <- fread(file, header=T, quote="\"", sep=",") ngs$time_numeric <- as.numeric(fastPOSIXct(ngs$Time)) ngs$game_play_id <- paste(ngs$GameKey, ngs$PlayID, sep="_") ngs<-ngs[, `:=`('TimeFromSnap'=time_numeric-time_numeric[Event=='ball_snap'][1]), by="game_play_id"] # head(ngs) # hist(ngs$TimeFromSnap) #extract starting points of all players, and where the punt was kicked, and where punt landed ngs[which(ngs$Event%in%c("punt_received","fair_catch","punt_land","punt_downed","touchback", "punt", "tackle", "ball_snap")| ngs$TimeFromSnap==2 | ngs$TimeFromSnap==3),] } ngs<-rbindlist(lapply(list_of_files, extractNGS)) ngs.organized<- ngs.organized<-ngs.organized[order(ngs.organized$time_numeric, decreasing = F),] ngs.organized<-ngs.organized[ngs.organized$Event%in% c('punt_received', 'fair_catch'),] ngs.organized<-ngs.organized[!duplicated(ngs.organized[, c("GSISID", "game_play_id")]),] all_player_roles<-merge(all_player_roles[, !colnames(all_player_roles)%in% c("x", "y"),with=F], ngs.organized[, c("GSISID","x", "y", "game_play_id"), with=F], by=c("GSISID","game_play_id"), all.x=T) all_player_roles<-all_player_roles[, `:=`(ball.loc.y=y[Role=="PLS"], ball.loc.x=x[Role=="PLS"]), by=c("game_play_id")] all_player_roles<-all_player_roles[, `:=`(pr.loc.y=y[Role=="PR"][1], pr.loc.x=x[Role=="PR"][1]), by=c("game_play_id")] all_player_roles$dist.to.pr<-sqrt((all_player_roles$y-all_player_roles$pr.loc.y)^2+(all_player_roles$x-all_player_roles$pr.loc.x)^2) ``` I also look at maximum impact, which I defined as the maximum decrease in speed a player experiences over the course of 2 seconds during a play. This too is significant for concussed players: ```{r include=T, eval=T , echo=F, message=F, warning=F} ggplot(all_player_roles, aes(x = max.impact, fill = concussion)) + geom_density(alpha=0.5, aes(fill=factor(concussion))) + labs(title="Density of concussion, given Max-Impact") + scale_x_continuous(breaks = scales::pretty_breaks(n = 10)) + theme_grey() ``` ## 2. Estimating Safety Improvements from Different Rule Changes For my analysis, my goal was to either limit max.impact or limit returns, as both are correlated to concussions. I decided to focus on modeling returns, as that was easier to predict. ```{r include=T, eval=T , echo=F} fit<-glm(Returned~FB+I(jammers>=3)+I(jammers>=4)+Box.Defender.dist.at.punt.x.max, family="binomial", data=play_information[!is.na(play_information$Box.Defender.dist.at.punt.x.max) ,]) summary(fit) ``` The above model estimates the effect of different variables on a punt being returned. Jammers and PFB on the return team seem to be most relevant to whether a kick is returned. Another variable I used is players dropping back into coverage. Many times a linebacker on the return team will fall back to block the gunner instead of rushing to block the punt. To account for this, I calculate the maximum distance beyond the line of scrimmage that a non-jammer or non-PFB player is at the time of snap. Finally, using these fitted coefficients, I can see the marginal effect of applying these different rules on the 2016-17 data. ```{r include=T, eval=T , echo=F} test1<-play_information[!is.na(play_information$Box.Defender.dist.at.punt.x.max),] test1$Group<-"Current Return Rate" test1$predReturned<-predict(fit, newdata=test1, type="response") test2<-play_information[!is.na(play_information$Box.Defender.dist.at.punt.x.max),] test2$jammers<-sapply(test2$jammers, function(x) min(x, 3)) test2$predReturned<-predict(fit, newdata=test2, type="response") test2$Group<-"At most 3 Jammers" test3<-play_information[!is.na(play_information$Box.Defender.dist.at.punt.x.max),] test3$jammers<-sapply(test3$jammers, function(x) min(x, 2)) test3$predReturned<-predict(fit, newdata=test3, type="response") test3$Group<-"At most 2 Jammers" test4<-play_information[!is.na(play_information$Box.Defender.dist.at.punt.x.max),] test4$Box.Defender.dist.at.punt.x.max<-sapply(test4$Box.Defender.dist.at.punt.x.max, function(x) min(x, 5)) test4$FB<-0 test4$jammers<-sapply(test4$jammers, function(x) min(x, 2)) test4$predReturned<-predict(fit, newdata=test4, type="response") test4$Group<-"At most 2 Jammers, \n box.defenders within 5 Yds of line, \n 0 PFB players" test<-do.call("rbind", list(test1, test2, test3, test4)) %>% data.table() test<-test[, list('Percent.Returned'=mean(predReturned)), by="Group"] test<-test[order(test$Percent.Returned),] ggplot(test, aes(x=Group, y=Percent.Returned)) + geom_bar(position=position_dodge(), stat="identity", colour="black", # Use black outlines, size=.3) ``` As you can see, all of these variables have an impact on a punt being returned. Switching from 3 to 2 jammers has the largest effect because of the frequency >=3 jammers occurs and because it has a large fitted coefficient. Finally, when implementing a rule, it's important to consider how teams will react. Number of PFB/jammers/or an LB dropping back are all inter-related, and if you limit jammers, teams would likely have a linebacker drop back when the punt is kicked. For that reason, it would make the most sense to implement all of restrictions together, as you would prevent teams from taking advantage. ## Conclusion In conclusion, my solution is to place limitations on return team formations. This includes:
1. At most 2 "jammers", and 2. No more than 3 players beyond 5 yards of the line of scrimmage at time of a punt (restricts the PFB position and restricts players dropping back on punts. Still allows for players to attempt to recover for fake punts) ** With my restrictions, the return rate would have gone down from 45% to 27%. This would have reduced the number of concussions by approximately:
1-(.27x6681x.010641950 +.73x6681x0.001592357)/(.45x6681x.010641950 +.55x6681x0.001592357)= 29%

Alternate solutions which would be worth analyzing more in depth are tighter restrictions on blindside blocking and head-to-head blocking, as well as changing the touchback line to the 25 yard line.