A Natural Gas Problem

Answer the questions to the best of your ability and be systematic in solving the problems. This will always take you far. An Excel Workbook called NaturalGas2019.xls and an identical R data file NaturalGas2019.RData contain the necessary data. Each sheet will contain data in both stacked and unstacked forms. Unstacked places the data for Henry Hub and Tianjin in separate data columns. Stacked data combines a single column of data with a single column describing which of the two locations the observation belongs to.

load(url("https://github.com/robertwwalker/DADMStuff/raw/master/NaturalGas2019.RData"))
library(tidyverse)
library(skimr)
library(ggplot2)

1. The central probability in a standard normal distribution between XXX and XXX is 0.6827.

result <- prob_norm(mean = 0, stdev = 1, plb = 0.159, pub = 0.841)
summary(result, type = "probs")

Probability calculator
Distribution: Normal
Mean        : 0 
St. dev     : 1 
Lower bound : 0.159 
Upper bound : 0.841 

P(X < -0.999) = 0.159
P(X > -0.999) = 0.841
P(X < 0.999) = 0.841
P(X > 0.999) = 0.159
P(-0.999 < X < 0.999)     = 0.682
1 - P(-0.999 < X < 0.999) = 0.318

plot(result, type = "probs")

0.6827 is the empirical rule probability between -1 and 1.

!Figure](https://github.com/robertwwalker/DADMStuff/raw/master/NormQ1.png)

result <- prob_norm(mean = 0, stdev = 1, plb = 0.1)
summary(result, type = "probs")

Probability calculator
Distribution: Normal
Mean        : 0 
St. dev     : 1 
Lower bound : 0.1 
Upper bound : 1 

P(X < -1.282) = 0.1
P(X > -1.282) = 0.9

plot(result, type = "probs")

The probability of a standard normal above -1.282 is 0.9.

3. For a normal distribution with mean 10 and standard deviation 5, what is the value of this normal such that 80% of the probability lies above it?

result <- prob_norm(mean = 10, stdev = 5, plb = 0.2)
summary(result, type = "probs")

Probability calculator
Distribution: Normal
Mean        : 10 
St. dev     : 5 
Lower bound : 0.2 
Upper bound : 1 

P(X < 5.792) = 0.2
P(X > 5.792) = 0.8

plot(result, type = "probs")

80 percent of the probability lies above 5.792.

The Case

A chemical manufacturer will build a new epoxy plant in one of two locations – one adjacent to an extant facility in Tianjin, China and one near the company’s headquarters in the United States near Henry Hub in Louisiana. The plant is a first-of-its-kind entirely automated plant that hosts only a tiny number of supervisors and labor costs across the two locations do not differ. The only input that might reasonably differ in cost is energy. Where to locate the factory and potential differences in the cost of energy will play a crucial role in our decision.

An Energy Cost Comparison

An intern has collected two random samples of energy prices from local newspapers for comparison. Over the last six months, the intern has a random sample of size 25 from Tianjin and a random sample of size 36 from Henry Hub reported in Newspapers – a sheet in the Excel Workbook on WISE referenced above. You are to undertake an initial summary of these data. You should note that both prices are reported in United States dollars per Million BTU (million British thermal units). The dates that the prices (spot) were recorded is listed in the adjacent column for each series. The worksheet contains the data in both a stacked and an unstacked format.

1. Provide a single plot combining the two samples of data. Which of the two series is more symmetric? What information in the plot led you to this choice?

visualize(
  Newspapers, 
  xvar = "Location", 
  yvar = "Price", 
  type = "box", 
  custom = FALSE
)

visualize(
  Newspapers, 
  xvar = "Price", 
  type = "density", 
  fill = "Location", 
  custom = FALSE
)

Henry Hub is both more variable and more symmetric. Tianjin has a fraction of the variability but is quite skewed.

2. What is the sample mean and sample standard deviation of natural gas prices in Tianjin? Henry Hub?

result <- explore(
  Newspapers, 
  vars = "Price", 
  byvar = "Location", 
  fun = c(
    "n_obs", "mean", "sd", "se", "min", "p25", "median", "p75", 
    "max"
  ), 
  nr = 2
)
summary(result)

Explore
Data        : Newspapers 
Grouped by  : Location 
Functions   : n_obs, mean, sd, se, min, p25, median, p75, max 
Top         : Function 

 Location variable n_obs  mean    sd    se   min   p25 median   p75   max
 HenryHub    Price    36 3.796 0.140 0.023 3.520 3.690  3.795 3.902 4.070
  Tianjin    Price    25 3.994 0.051 0.010 3.880 3.980  4.000 4.040 4.050

# dtab(result) %>% render()

The sample mean at Tianjin is $3.9944 with a standard deviation of $0.05 and a standard error of the mean of $0.01. The mean price at Henry Hub is $3.796 with a standard deviation of $0.14 and and standard error of the mean price of $0.023.

3. Assume that the underlying population is normal with the population mean and population standard deviation exactly equal to the sample mean and standard deviation that you have estimated from the sample from Tianjin. What is the probability of a daily price of $4.00 or higher at Tianjin?

result <- prob_norm(mean = 3.9944, stdev = 0.051, ub = 4)
summary(result)

Probability calculator
Distribution: Normal
Mean        : 3.994 
St. dev     : 0.051 
Lower bound : -Inf 
Upper bound : 4 

P(X < 4) = 0.544
P(X > 4) = 0.456

plot(result)

Assuming prices at Tianjin are normal with mean 3.9944 and standard deviation 0.051, the probability of a daily price of $4.00 or greater is 0.456.

4. What is a 95% confidence interval for the average price in Tianjin given this precise sample? Provide both boundaries.

With 95% confidence, the average price at Tianjin ranges from $3.97
result <- single_mean(Newspapers, var = "TianjinPrice")
summary(result)
plot(result, plots = "hist", custom = FALSE)

Error: <text>:1:6: unexpected numeric constant
1: With 95
         ^

With 95% confidence, the average price at Tianjin ranges from $3.973 to $4.015 per mBTU.

5. What is the probability of a mean price of $4.00 or lower given the Tianjin data?

result <- single_mean(
  Newspapers, 
  var = "TianjinPrice", 
  comp_value = 4, 
  alternative = "greater"
)
summary(result)

Single mean test
Data      : Newspapers 
Variable  : TianjinPrice 
Confidence: 0.95 
Null hyp. : the mean of TianjinPrice = 4 
Alt. hyp. : the mean of TianjinPrice is > 4 

  mean    sd n_obs n_missing
 3.994 0.051    25        36

   diff   se t.value p.value df    5% 100%  
 -0.006 0.01  -0.553   0.707 24 3.977  Inf  

Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

plot(result, plots = "hist", custom = FALSE)

The probability of a mean price of $4.00 or less is 0.707. It has to be greater than 0.5 because the sample mean is less than $4.00.

6. Speculative: What is a 95% confidence interval of the difference in average prices.

result <- compare_means(
  Newspapers, 
  var1 = "Location", 
  var2 = "Price"
)
summary(result, show = TRUE)

Pairwise mean comparisons (t-test)
Data      : Newspapers 
Variables : Location, Price 
Samples   : independent 
Confidence: 0.95 
Adjustment: None 

 Location  mean  n    sd    se    me
 HenryHub 3.796 36 0.140 0.023 0.047
  Tianjin 3.994 25 0.051 0.010 0.021

 Null hyp.            Alt. hyp.                       diff   p.value se    t.value df     2.5%  97.5%     
 HenryHub = Tianjin   HenryHub not equal to Tianjin   -0.199 < .001  0.025 -7.814  47.036 -0.25 -0.147 ***

Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

plot(result, plots = "scatter", custom = FALSE)

There is an estimated difference in the average prices that ranges between $0.147 and $0.25 per mBTU. Henry Hub is almost surely cheaper.

Profit Forecasting by Location FOR NEXT TIME

The division that is principally entrusted with forecasting and the supply chain has calculated a series of profit forecasts in dollars for the first 52 weeks after completion of the factory in Tianjin and the factory in Henry Hub. Because the production process is one-of-a-kind, it is known that the first 52 weeks will involve planned calibration of the factory; this is accounted for in the weekly forecasts and is the same pattern across both locations. These forecasts appear in worksheet Forecast of the workbook on WISE.

1. Summarize profits by location numerically and graphically. Try as best you can to capture the essence of the comparison in a single informative graphic.

The key here is capturing the relationship between the forecast periods and the forecasts themselves. There is a decided pattern. Let’s assume that we ignore that pattern. Let’s look at the boxplots.

visualize(
  Forecast, 
  xvar = "Location", 
  yvar = "Forecast", 
  type = "box", 
  custom = FALSE
)

The two distributions mostly overlap though one appears to be a slight shift of the other; Henry Hub being larger. But hidden beneath that is the fact that Henry Hub has forecast values that always exceed those for Tianjin. Let’s show that.

visualize(
  Forecast, 
  xvar = "Period", 
  yvar = "Forecast", 
  type = "line", 
  color = "Location", 
  custom = FALSE
)

result <- explore(
  Forecast, 
  vars = "Forecast", 
  byvar = "Location", 
  fun = c(
    "n_obs", "mean", "sd", "se", "min", "p25", "median", "p75", 
    "max"
  ), 
  nr = 2
)
summary(result)

Explore
Data        : Forecast 
Grouped by  : Location 
Functions   : n_obs, mean, sd, se, min, p25, median, p75, max 
Top         : Function 

 Location variable n_obs          mean            sd          se       min           p25        median           p75       max
 HenryHub Forecast    52 4,204,474.865 1,590,364.990 220,543.943 1,712,033 2,726,455.750 4,101,823.500 5,520,862.250 7,166,874
  Tianjin Forecast    52 3,860,590.442 1,561,618.097 216,557.466 1,525,784 2,336,259.000 3,713,691.000 5,197,883.000 6,761,850

# dtab(result) %>% render()

Henry Hub has average forecast profits of 4.2 million while Tianjin has an average forecast profit level of 3.86 million. The standard deviations and standard errors of the average are similar.

2. If profits in Henry Hub are normally distributed with mean and standard deviation as calculated above, what is the probability of profits between 4 million and 4.5 million?

result <- prob_norm(
  mean = 4204474.865, 
  stdev = 1590365, 
  lb = 4000000, 
  ub = 4500000
)
summary(result)

Probability calculator
Distribution: Normal
Mean        : 4204475 
St. dev     : 1590365 
Lower bound : 4000000 
Upper bound : 4500000 

P(X < 4000000) = 0.449
P(X > 4000000) = 0.551
P(X < 4500000) = 0.574
P(X > 4500000) = 0.426
P(4000000 < X < 4500000)     = 0.125
1 - P(4000000 < X < 4500000) = 0.875

plot(result)

Assuming a normal with exactly the sample mean and standard deviation that we estimated, we find a probability of 0.125 of being between 4 and 4.5 million dollars.

3. If profits in Tianjin are normally distributed with mean and standard deviation as calculated above, what is the range of profits to be obtained with central probability 0.5 [the interquartile range]?

result <- prob_norm(
  mean = 3860590.442, 
  stdev = 1561618.097, 
  plb = 0.25, 
  pub = 0.75
)
summary(result, type = "probs")

Probability calculator
Distribution: Normal
Mean        : 3860590 
St. dev     : 1561618 
Lower bound : 0.25 
Upper bound : 0.75 

P(X < 2807295.042) = 0.25
P(X > 2807295.042) = 0.75
P(X < 4913885.842) = 0.75
P(X > 4913885.842) = 0.25
P(2807295.042 < X < 4913885.842)     = 0.5
1 - P(2807295.042 < X < 4913885.842) = 0.5

plot(result, type = "probs")

Assuming a normal with exactly this mean and standard deviation, the middle 50% of forecast values ranges from $2.8 million and $4.9 million per period in forecast profits.

4. Provide a two-sided 90% confidence interval for average profits per week over the 52 week period in Tianjin, China.

result <- single_mean(
  Forecast, 
  var = "TianjinForecast", 
  comp_value = 3000000, 
  conf_lev = 0.9
)
summary(result)

Single mean test
Data      : Forecast 
Variable  : TianjinForecast 
Confidence: 0.9 
Null hyp. : the mean of TianjinForecast = 3000000 
Alt. hyp. : the mean of TianjinForecast is not equal to 3000000 

          mean            sd n_obs n_missing
 3,860,590.442 1,561,618.097    52        52

     diff       se t.value p.value df      5%     95%    
 860590.4 216557.5   3.974  < .001 51 3497795 4223386 ***

Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

plot(result, plots = "hist", custom = FALSE)

Average forecast profits range, with 90% confidence, between $3.5 million and $4.22 million per period for Tianjin.

5. Provide a two-sided 90% confidence interval for average profits per week over the 52 week period in Henry Hub.

result <- single_mean(
  Forecast, 
  var = "HenryHubForecast", 
  comp_value = 4, 
  conf_lev = 0.9
)
summary(result)

Single mean test
Data      : Forecast 
Variable  : HenryHubForecast 
Confidence: 0.9 
Null hyp. : the mean of HenryHubForecast = 4 
Alt. hyp. : the mean of HenryHubForecast is not equal to 4 

          mean            sd n_obs n_missing
 4,204,474.865 1,590,364.990    52        52

    diff       se t.value p.value df      5%     95%    
 4204471 220543.9  19.064  < .001 51 3835001 4573949 ***

Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

plot(result, plots = "hist", custom = FALSE)

Average forecast profits range from $3.835 million to $4.574 million per period for Henry Hub. It is worth noting the substantial overlap in likely values of the average. This is at odds with the visual evidence but it also shows how much the variability of the mean can be impacted by the volatility of the forecast pattern. These two forecasts are not and cannot be treated as independent samples.

6. Speculative: Provide a two-sided 95% confidence interval for the average difference in profits per week over the 52 week period..

It is an average difference; it is paired. The standard error of the difference is $13854 but the average difference is $343884. The t says there are 24.8 standard errors between the average difference and zero; there is a clear difference between the two that favors Henry Hub. Henry Hub is forecast to have $316072 to $371697 higher profits per period with 95% confidence.

result <- compare_means(
  Forecast, 
  var1 = "Location", 
  var2 = "Forecast", 
  samples = "paired"
)
summary(result, show = TRUE)

Pairwise mean comparisons (t-test)
Data      : Forecast 
Variables : Location, Forecast 
Samples   : paired 
Confidence: 0.95 
Adjustment: None 

 Location          mean  n            sd          se          me
 HenryHub 4,204,474.865 52 1,590,364.990 220,543.943 442,760.441
  Tianjin 3,860,590.442 52 1,561,618.097 216,557.466 434,757.255

 Null hyp.            Alt. hyp.                       diff     p.value se       t.value df 2.5%     97.5%       
 HenryHub = Tianjin   HenryHub not equal to Tianjin   343884.4 < .001  13853.59 24.823  51 316072.2 371696.7 ***

Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

plot(result, plots = "scatter", custom = FALSE)

7. Speculative: Is either of the two locations forecast to have higher profits with at least 95% confidence? Which one?

Henry Hub.

8. A bonus: If the building cost were $12 million more for the more profitable location, provide a 95% confidence interval for average recovery time of cost difference expressed in weeks.

How much do we retire of the 12 million dollar difference per week? 316072.2 to 371696.7 per week. So it takes 12 million divided by 316072.2 in the worst case and 12 million divided by 371696.7 per week in the best case yields and 32 to 38 week period to retire the difference.

Expert Assessment of Location

As management closes on a decision, word has leaked about the factory and the location decision. Boomblerg Television and BCNC have reported on the plans and industry experts now well understand the decision facing management. The Sall Tweet Kernel conducted a survey of independent samples of industry experts regarding the decision and their responses are recorded as worksheet Experts in the workbook on WISE. A Build represents a belief that the facility should be built in the given location (Tianjin/Henry Hub) while a No represents an opinion that the facility should not be built in the given location. With reference to these data, you are to answer the following questions.

1. Tablulate opinions by location. What is the proportion of analysts favoring each location? Provide some visualization relevant to at least one of the following questions and explain why.

result <- cross_tabs(
  Experts, 
  var1 = "Location", 
  var2 = "Build"
)
summary(result, check = c("observed", "row_perc"))

Cross-tabs
Data     : Experts 
Variables: Location, Build 
Null hyp.: there is no association between Location and Build 
Alt. hyp.: there is an association between Location and Build 

Observed:
         Build No  Total
HenryHub  66    15  81  
Tianjin   31    33  64  
Total     97    48 145  

Row percentages:
         Build   No Total
HenryHub  0.81 0.19     1
Tianjin   0.48 0.52     1
Total     0.67 0.33     1

Chi-squared: 17.628 df(1), p.value < .001 

0.0 % of cells have expected values below 5

plot(result, check = c("observed", "row_perc"), custom = FALSE)

2. Do analysts think the facility should be built at Henry Hub?

Yes. 81 percent do.

3. How large a sample would be required, with 95% confidence, to estimate the true proportion to within 0.1 assuming, for planning purposes, that the true proportion is 0.5?

result <- sample_size(type = "proportion")
summary(result)

Sample size calculation
Calculation type     : Proportion
Acceptable Error     : 0.1 
Proportion           : 0.5 
Confidence level     : 0.95 
Incidence rate       : 1 
Response rate        : 1 
Population correction: None

Required sample size     : 97
Required contact attempts: 97

4. The survey was, in both cases, sent to samples of size 100. As is clear, not everyone that was sent the survey actually responded. First, suppose that the probability of responding is 0.8. 95% of the time, we should receive XXX or more responses to this survey.

result <- prob_binom(n = 100, p = 0.8, pub = 0.05)
summary(result, type = "probs")

Probability calculator
Distribution: Binomial
n           : 100 
p           : 0.8 
Mean        : 80 
St. dev     : 4 
Lower bound : 
Upper bound : 0.05 (73)

P(X  = 73) = 0.022
P(X  < 73) = 0.034
P(X <= 73) = 0.056
P(X  > 73) = 0.944
P(X >= 73) = 0.966

plot(result, type = "probs")

5. Evaluate Henry Hub’s and Tianjin’s rates of response that you see in light of the above.

We have too little data on Tianjin with a response probability of 0.8.

6. Do analysts prefer either of the locations with any confidence?

result <- compare_props(
  Experts, 
  var1 = "Location", 
  var2 = "Build", 
  levs = "Build"
)
summary(result, show = TRUE)

Pairwise proportion comparisons
Data      : Experts 
Variables : Location, Build 
Level     : Build in Build 
Confidence: 0.95 
Adjustment: None 

 Location Build No  n     p    se    me
 HenryHub    66 15 81 0.815 0.043 0.085
  Tianjin    31 33 64 0.484 0.062 0.122

 Null hyp.            Alt. hyp.                       diff p.value chisq.value df 2.5%  97.5%    
 HenryHub = Tianjin   HenryHub not equal to Tianjin   0.33 < .001  17.628      1  0.182 0.479 ***

Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

plot(result, plots = "bar", custom = FALSE)

Yes. The probability of Build for Henry Hub is 0.182 to 0.479 higher than it is for Tianjin. Notice in the result, we see that the probability for Henry Hub is 0.815 with a margin of error of 0.085 while Tianjin is 0.484 with a margin of error of 0.122. So Henry Hub can be no lower than 0.73 and Tianjin is, at most, 0.606. There is an obvious difference. To quantify it, we turn to the confidence interval for the difference as described above.

Weather Concerns

Both locations [Tianjin and Henry Hub] are adjacent to large bodies of water. If the probability of more than 3 hurricanes in a year exceeds 0.05 in either location, this is a game changing fact. Management is unwilling to build in a location with such a high probability of a disatrous loss of investment and possible collateral damage.

1. If the arrival rate of hurricanes at Henry Hub is 1 per year (as it is known to be), what is the probability of more than 3 hurricanes in a given year at Henry Hub?

result <- prob_pois(lambda = 1, ub = 3)
summary(result)

Probability calculator
Distribution: Poisson
Lambda      : 1 
Mean        : 1 
Variance    : 1 
Lower bound :  
Upper bound : 3 

P(X  = 3) = 0.061
P(X  < 3) = 0.92
P(X <= 3) = 0.981
P(X  > 3) = 0.019
P(X >= 3) = 0.08

plot(result)

The result is less than 0.05.

2. If the arrival rate of hurricanes at Tianjin is 0.6 (as it is known 0.6o be), then what is the probability of more than 3 hurricanes in a given year at Tianjin?

result <- prob_pois(lambda = 0.6, ub = 3)
summary(result)

Probability calculator
Distribution: Poisson
Lambda      : 0.6 
Mean        : 0.6 
Variance    : 0.6 
Lower bound :  
Upper bound : 3 

P(X  = 3) = 0.02
P(X  < 3) = 0.977
P(X <= 3) = 0.997
P(X  > 3) = 0.003
P(X >= 3) = 0.023

plot(result)

The result is less than 0.05.

Final Problem

Defend a decision to build the factory in one of the two locations and relate each previous section and your analysis of the relevant data to your chosen alternative. If necessary, weigh the value of conflicting evidence in light of the source. You are only to faithfully report the results of your analysis and are to bring no external considerations to bear beyond the information in the data and the issues raised in this investigation.

Henry Hub Every bit of data supports that decision except the fact that hurricanes are more likely. But not that much so. The payback period is less than a year and that difference should yield benefits through the life of the project that exceed 300,000 dollars per week with 95% confidence. Natural gas is cheaper fifteen to twenty five cents cheaper, on average, and the experts think it is a better idea. There is no data to oppose building at Henry Hub.