% Copyright 2005 William S. Krasker
% W gives the probability (of scoring or winning) at the moment the whistle blows, if clock is running.
% S gives the probability (of scoring or winning) at the moment the ball is snapped for a
% real play (not a spike). This is also the win probability when the whistle blows if the clock stops.
% scorefactor=0 if offense trails by 4 points or more. The outputs of the model are then
% the probability of scoring a touchdown.
% scorefactor=1 if offense trails by 2 points or less. The outputs of the model are then
% the probability of winning.
% scorefactor=.5 if offense trails by 3 points. The outputs of the model are then
% the probability of winning.
scorefactor=1;
options = optimset('Display','off');
% probability of making a field goal from each field position.
fgprob=[.99 .98 .93 .85 .8 .74 .63 .5 .05 .01 zeros(1,11)];
spiketime=10; % Seconds required, after whistle blows, to spike the ball.
timeouttime=2; % Seconds required, after whistle blows, to call timeout.
letclockruntime=16; % Seconds required, after whistle blows, to run another real play.
probob=.2; % Probability that runner gets out of bounds after a completed pass.
maxgain=10; % This is in "units." Gain of 1 unit is 5 yards.
psacklose0= .20; %Probability that a sack loses 0 yards.
psacklose5= .39; %Probability that a sack loses 5 yards.
psacklose10=.41; %Probability that a sack loses 10 yards.
pstar=zeros(1,maxgain);
pstar(1)=.88;
pstar(5)=.35;
theta=1.6;
gamma=(pstar(1)/pstar(5)-1)/4^theta;
for gain=1:maxgain
pstar(gain)=pstar(1)/(1+gamma*(gain-1)^theta);
end
alpha=zeros(1,maxgain);
alpha(1)=.16;
alpha(maxgain)=.42;
for gain=2:maxgain-1
alpha(gain)=alpha(1)+((alpha(maxgain)-alpha(1))/(maxgain-1))*(gain-1);;
end
beta=.11; % Probability of an interception divided by the probability of an incompletion.
totalcoverage=7.5; % The constraint variable X.
psack=.0683; % Probability of a sack.
playtime=[7 7 7 7 7 7 7 7 7 7]; % Seconds used by a play, up until the whistle blows.
% Boundary conditions. If offense reaches goal line, they score.
S=zeros(21,121,5,4,4);
W=zeros(21,121,5,4,4);
for timeouts=1:4
for down=1:5
for togo=1:4
for t=1:121
S(1,t,down,togo,timeouts)=1;
W(1,t,down,togo,timeouts)=1;
end
end
end
end
% Backward induction
for t=2:121
for yards=2:20
for down=1:4
for togo=1:4
for timeouts=1:4
numgains=min(maxgain,yards-1);
sack=psacklose0*W(yards, max(t-playtime(1),1), down+1, togo, timeouts) ...
+psacklose5*W(min(yards+1,21), max(t-playtime(1),1), down+1, min(togo+1,4), timeouts) ...
+psacklose10*W(min(yards+2,21), max(t-playtime(1),1), down+1, min(togo+2,4), timeouts);
for gain=1:numgains
I(gain)=S(yards,max(t-playtime(gain),1),down+1,togo,timeouts);
Ibeta=I(gain)/(1+beta);
if gain >= togo
C(gain)=probob*S(yards-gain,max(t-playtime(gain),1),1,2,timeouts) ...
+(1-probob)*W(yards-gain,max(t-playtime(gain),1),1,2,timeouts);
else
C(gain)=probob*S(yards-gain,max(t-playtime(gain),1),down+1,togo-gain,timeouts) ...
+(1-probob)*W(yards-gain,max(t-playtime(gain),1),down+1,togo-gain,timeouts);
end
A1(gain)=pstar(gain)*(C(gain)-Ibeta)+Ibeta;
A2(gain)=alpha(gain)*Ibeta;
pbar(gain)=pstar(gain)/(1+alpha(gain)*(totalcoverage+0.1));
A3(gain)=pbar(gain)*(C(gain)-Ibeta)+Ibeta;
end
v(1)=max( C(1:numgains) );
v(2)=max( A3(1:numgains) );
if v(1)==0
S(yards,t,down,togo,timeouts)=fgprob(yards)*scorefactor;
else
[val,fvalfromfzero,exitflag]=fzero('sumfunc',v,options,numgains,A1,A2,alpha,totalcoverage);
if exitflag<0
exitflag
end
S(yards,t,down,togo,timeouts)=max( psack*sack+(1-psack)*val , fgprob(yards)*scorefactor );
end
spike=S(yards,max(t-spiketime,1),down+1,togo,timeouts);
if timeouts>1
calltime=S(yards,max(t-timeouttime,1),down,togo,timeouts-1);
else
calltime=0;
end
letclockrun=S(yards,max(t-letclockruntime,1),down,togo,timeouts);
W(yards,t,down,togo,timeouts)=max([spike calltime letclockrun]);
end
end
end
end
t
end