Page 8, Question 10

Your grandparents have an annuity. The value of the annuity increases each month by an automatic deposit of 1% interest on the previous month’s balance. Your grandparents withdraw $1000 at the beginning of each month for living expenses. Currently, they have $50,000 in the annuity. Model the annuity with a dynamical system.

\[ a_{n+1} = 1.01a_n - 1000 \]

Will the annuity run out of money? When?

# N0
amount = 50000
# Counter for the loop
i <- 0

# loop through periods in  dynamical model while balance > 0 
while (amount > 0){
    amount = 1.01*amount - 1000
    print(amount)
    i <- i + 1 
    print(i)
    if (amount < 1000) break # break out of loop if balance is less than 1000
}
## [1] 49500
## [1] 1
## [1] 48995
## [1] 2
## [1] 48484.95
## [1] 3
## [1] 47969.8
## [1] 4
## [1] 47449.5
## [1] 5
## [1] 46923.99
## [1] 6
## [1] 46393.23
## [1] 7
## [1] 45857.16
## [1] 8
## [1] 45315.74
## [1] 9
## [1] 44768.89
## [1] 10
## [1] 44216.58
## [1] 11
## [1] 43658.75
## [1] 12
## [1] 43095.34
## [1] 13
## [1] 42526.29
## [1] 14
## [1] 41951.55
## [1] 15
## [1] 41371.07
## [1] 16
## [1] 40784.78
## [1] 17
## [1] 40192.63
## [1] 18
## [1] 39594.55
## [1] 19
## [1] 38990.5
## [1] 20
## [1] 38380.4
## [1] 21
## [1] 37764.21
## [1] 22
## [1] 37141.85
## [1] 23
## [1] 36513.27
## [1] 24
## [1] 35878.4
## [1] 25
## [1] 35237.18
## [1] 26
## [1] 34589.56
## [1] 27
## [1] 33935.45
## [1] 28
## [1] 33274.81
## [1] 29
## [1] 32607.55
## [1] 30
## [1] 31933.63
## [1] 31
## [1] 31252.97
## [1] 32
## [1] 30565.5
## [1] 33
## [1] 29871.15
## [1] 34
## [1] 29169.86
## [1] 35
## [1] 28461.56
## [1] 36
## [1] 27746.18
## [1] 37
## [1] 27023.64
## [1] 38
## [1] 26293.87
## [1] 39
## [1] 25556.81
## [1] 40
## [1] 24812.38
## [1] 41
## [1] 24060.51
## [1] 42
## [1] 23301.11
## [1] 43
## [1] 22534.12
## [1] 44
## [1] 21759.46
## [1] 45
## [1] 20977.06
## [1] 46
## [1] 20186.83
## [1] 47
## [1] 19388.7
## [1] 48
## [1] 18582.58
## [1] 49
## [1] 17768.41
## [1] 50
## [1] 16946.09
## [1] 51
## [1] 16115.55
## [1] 52
## [1] 15276.71
## [1] 53
## [1] 14429.48
## [1] 54
## [1] 13573.77
## [1] 55
## [1] 12709.51
## [1] 56
## [1] 11836.6
## [1] 57
## [1] 10954.97
## [1] 58
## [1] 10064.52
## [1] 59
## [1] 9165.165
## [1] 60
## [1] 8256.817
## [1] 61
## [1] 7339.385
## [1] 62
## [1] 6412.779
## [1] 63
## [1] 5476.907
## [1] 64
## [1] 4531.676
## [1] 65
## [1] 3576.992
## [1] 66
## [1] 2612.762
## [1] 67
## [1] 1638.89
## [1] 68
## [1] 655.2788
## [1] 69

Starting with month 1, the annuity will end up with $655.2788 in month 69.

Page 17, Question 9

The data in the accompanying table show the speed \(n\) (in increments of 5 mph) of an automobile and the associated distance \(a_n\) in feet required to stop it once the brakes are applied. For instance, \(n = 6\) (representing \(6\times 5 = 30mph\)) requires a stopping distances of \(a_6 = 47ft\)

  1. Calculate and plot the change \(\Delta a_n\) versus \(n\). Does the graph reasonably approximate a linear relationship?
library(ggplot2)

# DaTA frame for the data in table
bDistancedf <- data.frame(n=1:16, an=c(3,6,11,21,32,47,65,87,112,140,171,204,241,282,325,376))

# Calc Delta N
bDistancedf$deltaAn = c(NA,tail(bDistancedf$an,-1)-head(bDistancedf$an,-1))

#Plot
ggplot(bDistancedf, aes(x=n, y=deltaAn)) + geom_point()
## Warning: Removed 1 rows containing missing values (geom_point).

The graph does reasonably approximate a linear relationship.

  1. Based on your conclusions in part (a), find a difference equation model for the stopping distance data. Test your model by plotting the errors in the predicted values against \(n\). Discuss the appropriateness of the model.

Since this is linear, we’re trying to find a line. We can also assume that if \(n = 0\), the stopping distance would be \(0\), and therefore the line passes through the origin.

To approximate the slope of this line, lets average the slopes from the origin to each of the points. In other words, for each point, lets calculate \(\frac{\Delta a_n}{n}\)

# Calc Slope y/x
bDistancedf$slope <- bDistancedf$deltaAn / bDistancedf$n

slope <- mean(bDistancedf$slope, na.rm=TRUE)

We could easily find the least squares solution to this problem, forcing the intercept to be 0, but the textbook was solving these problems in an approximate way by calculating the slope

altSlope <- summary(lm(deltaAn ~ 0 + n, data=bDistancedf))$coefficients["n","Estimate"]

Lets plot the errors in the predicted values vs n

bDistancedf$estimate <- slope*bDistancedf$n

bDistancedf$error <- bDistancedf$estimate - bDistancedf$deltaAn

ggplot(bDistancedf, aes(x=n, y=error)) + geom_point()
## Warning: Removed 1 rows containing missing values (geom_point).

This graph suggests that there is a pattern to the errors. This model may not be appropriate. If we were just fitting \(\Delta a_n\) to \(n\) using a linear model including an intercept, we could perhaps end up with a better fit.

Page 34, Question 13

Consider the spreading of a rumor through a company of 1000 employees, all working in the same building. We assume that the spreading of a rumor is similar to the spreading of a contagious disease (example 3, section 1.2) in that the number of people hearing the rumor each day is proportional to the product of the number who have heard the rumor previously and the number who have not heard the rumor. This is given by:

\[ r_{n+1} = r_n + kr_n (1000 - n) \]

Where \(k\) is a parameter that depends on how fast the rumor spreads and \(n\) is the number of days. Assume \(k = 0.001\) and further assume that four people initially have heard the rumor. How soon will all 1000 employees have heard the rumor?

rumordf <- data.frame(n=0, r=4, roundr=4)

k <- 0.001

for(i in 1:13){
  rumordf <- rbind(rumordf, c(i, rumordf$r[i] + k*rumordf$r[i]*(1000 - rumordf$r[i]),
                              round(rumordf$roundr[i] + k*rumordf$roundr[i]*(1000 - 
                                                                               rumordf$roundr[i]))))
}

rumordf
##     n          r roundr
## 1   0    4.00000      4
## 2   1    7.98400      8
## 3   2   15.90426     16
## 4   3   31.55557     32
## 5   4   62.11538     63
## 6   5  120.37244    122
## 7   6  226.25535    229
## 8   7  401.31922    406
## 9   8  641.58132    647
## 10  9  871.53605    875
## 11 10  983.49701    984
## 12 11  999.72765   1000
## 13 12  999.99993   1000
## 14 13 1000.00000   1000

I created two columns here, one that rounds the number of people who have heard the rumor to the nearest integer. Given the initial condition of four people who have heard the rumor, this doesn’t make too big a difference. After 12 days, the rumor will have spread throughout the office.

Page 55, Question 6

An economist is interested in the variation of the price of a single product. It is observed that a high price for the product in the market attracts more suppliers. However, increasing the quantity of the product supplied tends to drive the price down. Over time, there is an interaction between price and supply. The economist has proposed the following model,p where \(P_n\) represents the price of the product at year \(n\), and \(Q_n\) represents the quantity. Find the equilibrium values for this system.

\[ P_{n+1} = P_n - 0.1(Q_n - 500) \] \[ Q_{n+1} = Q_n + 0.2(P_n - 100) \]

  1. Does the model make sense intuitively? What is the significance of the constants 100 and 500? Explain the significants of the constance -0.1 and 0.2.

See explanation below.

  1. Test the initial conditions in the following table and predict the long-term behavior:
Case Price Quantity
A 100 500
B 200 500
C 100 600
D 100 400
library(reshape2)

iterations <- 200

pinit <- 100
qinit <- 500

econdf <- data.frame(n = 0, p = pinit, q= qinit)

for(i in 1:iterations){
  econdf <- rbind(econdf,
                  c(i, tail(econdf$p,1) - 0.1*(tail(econdf$q,1) - 500),
                    tail(econdf$q,1) + 0.2*(tail(econdf$p,1) - 100)))
}

ggplot(melt(econdf, id="n"), aes(x=n, y=value, color=variable)) + geom_line()

#################################################################################

pinit <- 200
qinit <- 500

econdf <- data.frame(n = 0, p = pinit, q= qinit)

for(i in 1:iterations){
  econdf <- rbind(econdf,
                  c(i, tail(econdf$p,1) - 0.1*(tail(econdf$q,1) - 500),
                    tail(econdf$q,1) + 0.2*(tail(econdf$p,1) - 100)))
}

ggplot(melt(econdf, id="n"), aes(x=n, y=value, color=variable)) + geom_line()

#################################################################################

pinit <- 100
qinit <- 600

econdf <- data.frame(n = 0, p = pinit, q= qinit)

for(i in 1:iterations){
  econdf <- rbind(econdf,
                  c(i, tail(econdf$p,1) - 0.1*(tail(econdf$q,1) - 500),
                    tail(econdf$q,1) + 0.2*(tail(econdf$p,1) - 100)))
}

ggplot(melt(econdf, id="n"), aes(x=n, y=value, color=variable)) + geom_line()

#################################################################################

pinit <- 100
qinit <- 400

econdf <- data.frame(n = 0, p = pinit, q= qinit)

for(i in 1:iterations){
  econdf <- rbind(econdf,
                  c(i, tail(econdf$p,1) - 0.1*(tail(econdf$q,1) - 500),
                    tail(econdf$q,1) + 0.2*(tail(econdf$p,1) - 100)))
}

ggplot(melt(econdf, id="n"), aes(x=n, y=value, color=variable)) + geom_line()

Seeing the graphs made me realize what was going on intuitively. At first I assumed this system would result in an equilibrium price and quantity. While this is possible given the right initial conditions (a price of 100 and a quantity of 500), any other values will cause the system to oscillate without bound around the equilibrium.

The difference in sign between 0.1 and 0.2 in the two equations is what causes the oscillation. If these were both positive or both negative, they would increase or decrease without bound most of the time (this depends on what the initial price and quantity are).

The one equilibrium case I found was if the equations are the same, the second term is being subtracted from the first, and the initial price and quantity are the same:

pinit <- 20
qinit <- 20

econdf <- data.frame(n = 0, p = pinit, q= qinit)

for(i in 1:iterations){
  econdf <- rbind(econdf,
                  c(i, tail(econdf$p,1) - 0.1*(tail(econdf$q,1) - 500),
                    tail(econdf$q,1) - 0.1*(tail(econdf$p,1) - 500)))
}

ggplot(melt(econdf, id="n"), aes(x=n, y=value, color=variable)) + geom_line()

Any differences from this specific set of cases will once again cause divergence. For example, in the above case if the initial price is not equal to the initial quantity:

pinit <- 21
qinit <- 20

econdf <- data.frame(n = 0, p = pinit, q= qinit)

for(i in 1:iterations){
  econdf <- rbind(econdf,
                  c(i, tail(econdf$p,1) - 0.1*(tail(econdf$q,1) - 500),
                    tail(econdf$q,1) - 0.1*(tail(econdf$p,1) - 500)))
}

ggplot(melt(econdf, id="n"), aes(x=n, y=value, color=variable)) + geom_line()

This is probably not a good model for the interaction between price and quantity in normal markets. However, I could imagine some economic processes that might function in this way. This sort of oscillating behavior reminded me of the “fool in the shower” metaphor for monetary policy, where central banks with imperfect information react to the economy in a delayed fashion, raising and lowering interest rates in a way that’s out of sync with the rest of the economy.