*** To read a text file with data in multiple columns and rows and
    want to plot specific columns from the file. First store all data
    into a matrix, here it's 5 columns and 282 rows. Then read "All"
    rows 1 -> 282 and specific columns (x, y): 

*** Figure 1:
<<Graphics`MultipleListPlot`

m = ReadList["P:\My Documents\FrankMasci\Easter2\Liquor\Residuals.txt", Table[Number,{5}], 282]
 
data1 = m[[All, {4, 5}]]

g = MultipleListPlot[data1,
  PlotJoined -> {False},
  SymbolShape -> {PlotSymbol[Box, 1, Filled -> True]},
  SymbolStyle -> {Hue[0.6]},
  PlotRange -> {{0.9, 1.1}, {0.9, 1.1}},
  AxesOrigin -> {1, 1},
  FrameLabel -> {StyleForm["X12-ARIMA residual: exp[log(obs)-log(mod)]",
                 FontFamily -> "Times", FontSize -> 11],
                 StyleForm["X12-D13 irregular",
                 FontFamily -> "Times",FontSize -> 11]},
  RotateLabel -> False,
  Frame -> True]

h = Show[g, {Graphics[{Hue[1.0], Line[{{0.9, 0.9}, {1.1, 1.1}}]}]}]
SetDirectory["P:\My Documents\FrankMasci\Easter2"]
Export["Figure1.gif", h, ImageSize -> Automatic, ImageResolution -> 150]

*** Figure 2:
<<Graphics`MultipleListPlot`

m = ReadList["P:\My Documents\FrankMasci\Easter2\SuperMarket\Residuals.txt", Table[Number,{5}], 281]
 
data1 = m[[All, {4, 5}]]

g = MultipleListPlot[data1,
  PlotJoined -> {False},
  SymbolShape -> {PlotSymbol[Box, 1, Filled -> True]},
  SymbolStyle -> {Hue[0.6]},
  PlotRange -> {{0.9, 1.1}, {0.9, 1.1}},
  AxesOrigin -> {1, 1},
  FrameLabel -> {StyleForm["X12-ARIMA residual: exp[log(obs)-log(mod)]",
                 FontFamily -> "Times", FontSize -> 11],
                 StyleForm["X12-D13 irregular",
                 FontFamily -> "Times",FontSize -> 11]},
  RotateLabel -> False,
  Frame -> True]

h = Show[g, {Graphics[{Hue[1.0], Line[{{0.9, 0.9}, {1.1, 1.1}}]}]}]
SetDirectory["P:\My Documents\FrankMasci\Easter2"]
Export["Figure2.gif", h, ImageSize -> Automatic, ImageResolution -> 150]

================================================================================

*** Figure 3:
<<Graphics`MultipleListPlot`

m = ReadList["P:\My Documents\FrankMasci\Easter2\Liquor\Residuals.txt", Table[Number,{5}], 282]

data1 = m[[All, {1, 4}]]
data2 = m[[All, {1, 5}]]

g = MultipleListPlot[data1, data2,
  PlotJoined -> {True, True},
  PlotStyle -> {{Hue[0.6],Thickness[0.005]},
		{Hue[1.0],Thickness[0.005]}},
  SymbolShape -> {PlotSymbol[Box, 0.1, Filled -> True],
                  PlotSymbol[Box, 0.1, Filled -> True]},
  SymbolStyle -> {Hue[0.6], Hue[1.0]},
  PlotRange -> {{-5, 287}, {0.85, 1.15}},
  AxesOrigin -> {-5, 1},
  PlotLegend -> {StyleForm["Exp[ARIMA residual]",
                 FontFamily -> "Times", FontSize -> 11],
                 StyleForm["D13 Irregular",
                 FontFamily -> "Times",FontSize -> 11]},
  LegendBorder -> GrayLevel[1],
  LegendPosition -> {-1.3, 0.15},
  LegendSpacing -> -0.05,
  FrameLabel -> {StyleForm["Month# Mar'82 - Aug'05",
                 FontFamily -> "Times", FontSize -> 11],
                 StyleForm["Residual or irregular",
                 FontFamily -> "Times",FontSize -> 11]},
  RotateLabel -> False,
  Frame -> True]

SetDirectory["P:\My Documents\FrankMasci\Easter2"]
Export["Figure3.gif", g, ImageSize -> Automatic, ImageResolution -> 150]

================================================================================
*** histogram of raw unbinned input data - here it's column 3 of input file
    and normalize all values such that have PDF normalized to unity
    for comparison with known PDFs.

*** Figure 5:
<< Graphics`Graphics`

m = ReadList["P:\My Documents\FrankMasci\Easter2\Liquor\Residuals.txt", Table[Number,{5}], 282]

data1 = m[[All, 3]]

g = Histogram[data1,
              HistogramScale -> 1,
              HistogramRange-> {-0.11,0.1},
              Frame -> True,
              FrameLabel -> {StyleForm["Raw ARIMA residual",
	         FontFamily -> "Times", FontSize -> 11],
                 StyleForm["Number",
		 FontFamily -> "Times",FontSize -> 11]},
                 RotateLabel->False]

SetDirectory["P:\My Documents\FrankMasci\Easter2"]
Export["Figure5.gif", g, ImageSize -> Automatic, ImageResolution -> 150]

================================================================================
*** normalized histogram of raw data with overlayed pdf function.

<< Graphics`Graphics`
<< Statistics`NormalDistribution`
<< Statistics`DescriptiveStatistics`

m = ReadList["P:\My Documents\FrankMasci\Easter2\Liquor\Residuals.txt", Table[Number,{5}], 282]

data1 = m[[All, 3]]

g = Histogram[data1,
              HistogramScale -> 1,
              HistogramRange-> {-0.11,0.1},
              Frame -> True]

# sample stats:
v = Variance[data1]
Kurtosis[data1]
Skewness[data1]
LocationReport[data1]
DispersionReport[data1]
ShapeReport[data1]

n = 0.03;

d = NormalDistribution[0, v];

p = ListPlot[Table[{x,n*PDF[d, x]}, {x, -0.1, 0.1, 0.01}],
             PlotJoined->True] 

Show[g, p]

================================================================================

*** plot known pdfs/cdfs..
*** Under "ContinuousDistributions" in help browser"

*** NORMAL DISTRIBUTION: mu = 0, var = 1:
*** N.B. only for Gaussian: mu +/- 1s.d. contains 68.3% of values!

<<Statistics`NormalDistribution` 

d = NormalDistribution[0, 1];
med = Quantile[d, 0.5];
avg = Mean[d];
sd = StandardDeviation[d];
Pravg = N[PDF[d, avg]];
L = avg - sd;
U = avg + sd;
PrsdL = N[PDF[d, L]];
PrsdU = N[PDF[d, U]];

g = Plot[PDF[d, x], {x, -5, 5},
  AxesOrigin -> {-5, 0},
  FrameLabel -> {StyleForm["x", FontFamily -> "Times", FontSize -> 11],
                 StyleForm["P(x)", FontFamily -> "Times",FontSize -> 11]},
  RotateLabel -> False,
  Frame -> True]

Show[g, {Graphics[{Hue[0.6], Line[{{avg, 0}, {avg, Pravg}}]}],
         Graphics[{Hue[1.0], Line[{{L, 0}, {L, PrsdL}}]}],
         Graphics[{Hue[1.0], Line[{{U, 0}, {U, PrsdU}}]}]}]

*** Compute kurtosis of d. = 3 for Gaussian.
k = N[Kurtosis[d]]

*** Compute CDF (probability) within U, L values from above:
N[CDF[d, U]] - N[CDF[d, L]] = 0.682689

------------------------
*** LOGNORMAL DISTRIBUTION: mu = 0, var = 1 for normally distributed variable "log[x]".
*** does not have symmetric confidence intervals.

*** N.B must not have Statistics`NormalDistribution` simultaneously loaded.

<<Statistics`ContinuousDistributions`

d = LogNormalDistribution[0, 1];
med = Quantile[d, 0.5];
avg = Mean[d];
sd = StandardDeviation[d];
Pravg = N[PDF[d, avg]];
L = avg - sd;
U = avg + sd;
PrsdL = N[PDF[d, L]];
PrsdU = N[PDF[d, U]];

g = Plot[PDF[d, x], {x, 0, 5},
  AxesOrigin -> {-1, 0},
  FrameLabel -> {StyleForm["x", FontFamily -> "Times", FontSize -> 11],
                 StyleForm["P(x)", FontFamily -> "Times",FontSize -> 11]},
  RotateLabel -> False,
  Frame -> True]

Show[g, {Graphics[{Hue[0.6], Line[{{avg, 0}, {avg, Pravg}}]}],
         Graphics[{Hue[1.0], Line[{{L, 0}, {L, PrsdL}}]}],
         Graphics[{Hue[1.0], Line[{{U, 0}, {U, PrsdU}}]}]}]

*** Compute CDF (probability) within U, L (L = 0 here) values from above:
N[CDF[d, U]] - N[CDF[d, L]] = 0.909488


================================================================================

*** histogram of raw unbinned input data - here it's columns 4 and 5 of input file
    and normalize all values such that have PDF normalized to unity
    for comparison with known PDFs.

***Figure 6:
<< Graphics`Graphics`

m = ReadList["P:\My Documents\FrankMasci\Easter2\Liquor\Residuals.txt", Table[Number,{5}], 282]

data1 = m[[All, 4]]
data2 = m[[All, 5]]

g1 = Histogram[data1,
	      BarStyle -> {Hue[1.0]},
              HistogramScale -> 1,
              HistogramRange-> {0.9,1.1},
              Frame -> True]


g2 = Histogram[data2,
	      BarStyle -> {Hue[0.6]},
              HistogramScale -> 1,
              HistogramRange-> {0.9,1.1},
              Frame -> True]

h = Show[g1, g2,
  FrameLabel -> {StyleForm["Exp[residual] (red) or Irreg. (blue)",
	         FontFamily -> "Times", FontSize -> 11],
                 StyleForm["Number",
		 FontFamily -> "Times",FontSize -> 11]},
                 RotateLabel->False];

SetDirectory["P:\My Documents\FrankMasci\Easter2"]
Export["Figure6.gif", h, ImageSize -> Automatic, ImageResolution -> 150]

================================================================================

*** Matching with Logmormal distributions: i.e., above overlayed with best fitting
    lognormal distributions.

*** Figure 8:
<< Graphics`Graphics`
<< Statistics`ContinuousDistributions`

m = ReadList["P:\My Documents\FrankMasci\Easter2\Liquor\Residuals.txt", Table[Number,{5}], 282]

data1 = m[[All, 4]]
data2 = m[[All, 5]]

g1 = Histogram[data1,
	      BarStyle -> {Hue[1.0]},
              HistogramScale -> 1,
              HistogramRange-> {0.9,1.1},
              Frame -> True]


g2 = Histogram[data2,
	      BarStyle -> {Hue[0.6]},
              HistogramScale -> 1,
              HistogramRange-> {0.9,1.1},
              Frame -> True]

g = Show[g1, g2,
  FrameLabel -> {StyleForm["Exp[residual] (red) or Irreg. (blue)",
	         FontFamily -> "Times", FontSize -> 11],
                 StyleForm["Number",
		 FontFamily -> "Times",FontSize -> 11]},
                 RotateLabel->False];

d = LogNormalDistribution[0, 0.02];
p1 = ListPlot[Table[{x, 1.8*PDF[d, x]}, {x, 0, 3, 0.005}],
              PlotJoined -> True,
              PlotStyle -> {Hue[0.6], Thickness[0.01]}];
gn1 = Show[g, p1];

d = LogNormalDistribution[0, 0.025];
p2 = ListPlot[Table[{x, PDF[d, x]}, {x, 0, 3, 0.005}],
     PlotJoined -> True,
     PlotStyle -> {Hue[1.0], Thickness[0.01]}];
gn2 = Show[gn1, p2];

SetDirectory["P:\My Documents\FrankMasci\Easter2"]
Export["Figure8.gif", gn2, ImageSize -> Automatic, ImageResolution -> 150]

================================================================================

***Figure 7: comparing lognormal and normal pdfs for small variances.

<< Graphics`Legend`

v = 0.01;
frac1 = (((1/R)*Exp[-((Log[R])^2)/(2*v)])/(Exp[-((R - 1)^2)/(2*((Exp[v]-1)*Exp[v]))]));

v = 0.001;
frac2 = (((1/R)*Exp[-((Log[R])^2)/(2*v)])/(Exp[-((R - 1)^2)/(2*((Exp[v]-1)*Exp[v]))]));

v = 0.0001;
frac3 = (((1/R)*Exp[-((Log[R])^2)/(2*v)])/(Exp[-((R - 1)^2)/(2*((Exp[v]-1)*Exp[v]))]));

g1 = Plot[{frac1, frac2, frac3}, {R, 0.8,1.2},
         Frame -> True,
         TextStyle -> {FontSlant -> "Plain", FontSize -> 6},
         AxesOrigin -> {1, 1},
         PlotRange -> {{0.8, 1.2}, {0.8, 1.2}},
         PlotStyle -> {Hue[1.0],Hue[0.6],Hue[0.4]},
         FrameLabel -> {StyleForm["R",
                 FontFamily -> "Times", FontSize -> 10],
                 StyleForm["LogNormal(R)/Normal(R)",
                 FontFamily -> "Times",FontSize -> 10]},
         PlotLegend -> {StyleForm["var(R) = 0.01",
                        FontFamily -> "Times", FontSize -> 9],
                        StyleForm["var(R) = 0.001",
                        FontFamily -> "Times", FontSize -> 9],
                        StyleForm["var(R) = 0.0001",
                        FontFamily -> "Times", FontSize -> 9]},
         LegendBorder -> GrayLevel[1],
         LegendPosition -> {1.07,0.05},
         LegendSpacing -> -0.1,
	 LegendShadow -> None];

SetDirectory["P:\My Documents\FrankMasci\Easter2"]
Export["Figure7.gif", g1, ImageSize -> Automatic, ImageResolution -> 150]

================================================================================

Easter proximity charts using D13 irregulars (1 = neutral):

<<Graphics`MultipleListPlot`

m = ReadList["P:\My Documents\FrankMasci\Easter2\Liquor\EasterProx_D13.tbl", Table[Number,{4}], 23] 
data1 = m[[All, {2, 3}]]
data2 = m[[All, {2, 4}]] 

g = MultipleListPlot[data1, data2,
  PlotJoined -> {False, False},
  SymbolShape -> {PlotSymbol[Box, 3, Filled -> False],
                  PlotSymbol[Box, 3, Filled -> True]},
  SymbolStyle -> {Hue[0.6], Hue[1.0]},
  PlotRange -> {{-8, 25}, {0.94, 1.06}},
  AxesOrigin -> {0, 1},
  PlotLegend -> {StyleForm["MARCH", FontFamily -> "Times", FontSize -> 11],
                 StyleForm["APRIL", FontFamily -> "Times",FontSize -> 11]},
  LegendBorder -> GrayLevel[1],
  LegendPosition -> {1.07,0.05},
  LegendSpacing -> -0.3,
  FrameLabel -> {StyleForm["Easter date relative to Mar. 31", FontFamily -> "Times", FontSize -> 14],
                 StyleForm["D13 Irregular", FontFamily -> "Times",FontSize -> 14]},
  RotateLabel -> False,
  Frame -> True]

================================================================================

***Figure 9:

Easter proximity charts using **raw** ARIMA residuals (0 = neutral):

<<Graphics`MultipleListPlot`

m = ReadList["P:\My Documents\FrankMasci\Easter2\appliances\EasterProx_RSD.tbl", Table[Number,{4}], 23] 
data1 = m[[All, {2, 3}]]
data2 = m[[All, {2, 4}]] 

g = MultipleListPlot[data1, data2,
  PlotJoined -> {False, False},
  SymbolShape -> {PlotSymbol[Box, 3, Filled -> False],
                  PlotSymbol[Box, 3, Filled -> True]},
  SymbolStyle -> {Hue[0.6], Hue[1.0]},
  PlotRange -> {{-8, 25}, {-0.12, 0.12}},
  AxesOrigin -> {0, 0},
  PlotLegend -> {StyleForm["MARCH", FontFamily -> "Times", FontSize -> 11],
                 StyleForm["APRIL", FontFamily -> "Times",FontSize -> 11]},
  LegendBorder -> GrayLevel[1],
  LegendPosition -> {1.07,0.05},
  LegendSpacing -> -0.3,
  FrameLabel -> {StyleForm["Easter date relative to Mar. 31", FontFamily -> "Times", FontSize -> 14],
                 StyleForm["ARIMA Residual", FontFamily -> "Times",FontSize -> 14]},
  RotateLabel -> False,
  Frame -> True]

================================================================================

Simulation to explore what the Median Absolute Deviation (MAD) measure represents for
the Normal distribution. Answer: MAD = the 75 percentile for Normal distbn or C.I. of 50%
about the mean or median (since median = mean for Normal distbn.) Note: 1sigma for Normal
distribution = 68.3% C.I.

<< Statistics`ContinuousDistributions`
ndist = NormalDistribution[0, 1]
data = Table[Random[ndist], {5000}]
mad = Median[Abs[data - Median[data]]]
std = StandardDeviation[data]

Quantile[ndist, 0.75]
StandardDeviation[ndist]

================================================================================

***Figure 1: version 4 document:

Easter proximity charts using **raw** ARIMA residuals (0 = neutral):

<<Graphics`MultipleListPlot`

m = ReadList["P:\My Documents\FrankMasci\Easter2\Liquor_OrigwithSimOutliersWithEPcorr\EasterProx_RSD.tbl", Table[Number,{4}], 23] 
data1 = m[[All, {2, 3}]]
data2 = m[[All, {2, 4}]] 

g = MultipleListPlot[data1, data2,
  PlotJoined -> {False, False},
  SymbolShape -> {PlotSymbol[Box, 3, Filled -> False],
                  PlotSymbol[Box, 3, Filled -> True]},
  SymbolStyle -> {Hue[0.6], Hue[1.0]},
  PlotRange -> {{-8, 25}, {-0.2, 0.2}},
  AxesOrigin -> {0, 0},
  PlotLegend -> {StyleForm["MARCH", FontFamily -> "Times", FontSize -> 11],
                 StyleForm["APRIL", FontFamily -> "Times",FontSize -> 11]},
  LegendBorder -> GrayLevel[1],
  LegendPosition -> {1.07,0.05},
  LegendSpacing -> -0.3,
  FrameLabel -> {StyleForm["Easter date relative to Mar. 31", FontFamily -> "Times", FontSize -> 14],
                 StyleForm["ARIMA Residual", FontFamily -> "Times",FontSize -> 14]},
  RotateLabel -> False,
  Frame -> True]

================================================================================

***Figure 2: version 4 document:

Easter proximity charts using **raw** ARIMA residuals (0 = neutral):

<<Graphics`MultipleListPlot`

m = ReadList["P:\My Documents\FrankMasci\Easter2\Liquor_WithWithoutBestEPmod\WithBestEPmod\EasterProx_RSD.tbl", Table[Number,{4}], 23] 
data1 = m[[All, {2, 3}]]
data2 = m[[All, {2, 4}]] 

g = MultipleListPlot[data1, data2,
  PlotJoined -> {False, False},
  SymbolShape -> {PlotSymbol[Box, 3, Filled -> False],
                  PlotSymbol[Box, 3, Filled -> True]},
  SymbolStyle -> {Hue[0.6], Hue[1.0]},
  PlotRange -> {{-8, 25}, {-0.08, 0.08}},
  AxesOrigin -> {0, 0},
  PlotLegend -> {StyleForm["MARCH", FontFamily -> "Times", FontSize -> 11],
                 StyleForm["APRIL", FontFamily -> "Times",FontSize -> 11]},
  LegendBorder -> GrayLevel[1],
  LegendPosition -> {1.07,0.05},
  LegendSpacing -> -0.3,
  FrameLabel -> {StyleForm["Easter date relative to Mar. 31", FontFamily -> "Times", FontSize -> 14],
                 StyleForm["ARIMA Residual", FontFamily -> "Times",FontSize -> 14]},
  RotateLabel -> False,
  Frame -> True]

================================================================================