# 1. Read in lists of digit duration data.
#
# Each list is generated from the segmentation of 100 strings of 7 digits,
#    grouped in the style of an American phone number, 3 + 4,
#    e.g. "123-4567")

X4 <- read.table("MYL_NList4.out")
X8 <- read.table("MYL_List8.ldur1")

#   NList4 is arranged so that
#        each digit occurs 10 times in each position
#        and each pair of digits occurs once spanning each pair of positions
#   (List8 is not entirely balanced in that way)

# Lines are of the form
#
# #001 3 0.200 1 0.234 0 0.328 6 0.302 0 0.305 7 0.351 6 0.494
#
# where column 1 is the string number "#001" to "#100"
#   cols 2,4,6,8,10,12 are the digits in the string
#   and cols 3,5,7,9,11,13,15 are the durations
#
# So pull out the digits and durations for each table separately:


Nums4 <- X4[,c(2,4,6,8,10,12,14)]
Durs4 <- X4[,c(3,5,7,9,11,13,15)]

Nums8 <- X8[,c(2,4,6,8,10,12,14)]
Durs8 <- X8[,c(3,5,7,9,11,13,15)]

# 2. Plot overall mean durations, separately for each list of 100

png(filename="MYL_Digits%03d.png", width=480, height=480)

yr <- range(c(mean(Durs4), mean(Durs8)))

plot(1:7,mean(Durs4), type="b", pch="x", col="red", ylim=yr,
     xlab="Digit Position", ylab="Mean Duration (seconds)",
     main="Mean duration for all digits\n-x- NList4    -o- List8")

points(1:7,mean(Durs8), type="b", pch="o", col="blue")

# 3. Now plot the mean durations digit by digit
#
# First get the mean durations in a convenient form:

Dmeans <- matrix(nrow=10,ncol=7)

for (pos in 1:7){
  for(digit in 0:9){
    Dmeans[digit+1,pos] <- mean(Durs4[Nums4[,pos]==digit,pos])
  }
}

# ... and plot them:

cols <- c("red","darkgreen", "blue", "orange", "brown", "black", "maroon", 
           "purple", "rosybrown", "magenta" )
ltys <- c(1,2,3,4,5,6,1,2,3,4)
yr <- range(Dmeans)

plot(1:7, Dmeans[1,],
       type="n", ylim=yr,
       xlab="Digit Position", ylab="Mean Duration (seconds)",
       main="Mean Duration for Digits 0-9")

for (n in 0:9){
   points(1:7, Dmeans[n+1,],
       type="b", col=cols[n+1], lty=ltys[n+1], pch=as.character(n))
}

#
# 4. Now plot the relationship between duration in position 2 and position 3,
#     digit by digit:

lims <- range(c(Dmeans[,2], Dmeans[,3]))

plot(c(.3,.3),c(.3,.3), ylim=lims, xlim=lims, type="n",
     xlab="Position 2 Duration", ylab="Position 3 Duration",
     main="Internal Pre-boundary Lengthening, By Digit")

text(Dmeans[,2],Dmeans[,3], labels=0:9, col=cols)
abline(a=0,b=1, col="black", lty="dotted")

# Fit a line just to monosyllabic digits c(1,2,3,4,5,6,8,9) 

monosyl <- c(1,2,3,4,5,6,8,9)+1  # add 1 since digit 0 is row 1
p1 <- glm(Dmeans[monosyl,3] ~ Dmeans[monosyl,2])

# And plot it:
abline(p1[1], lty="dashed", col="green")

# 5. Ditto for position 6 and position 7:

lims <- range(c(Dmeans[,6], Dmeans[,7]))

plot(c(.3,.3),c(.3,.3), ylim=lims, xlim=lims, type="n",
     xlab="Position 6 Duration", ylab="Position 7 Duration",
     main="Final Pre-boundary Lengthening, By Digit")

text(Dmeans[,6],Dmeans[,7], labels=0:9, col=cols)
abline(a=0,b=1, col="black", lty="dotted")

# Add line for c(1,2,3,4,9)
#     i.e. those without the complications of 
#      - final unreleased stop (8)
#      - two syllables (0,7)
#      - extra lengthening due to final fricative (5,6)

which <- c(1,2,3,4,9)+1
p2 <- glm(Dmeans[which,7] ~ Dmeans[which,6])
abline(p2[1], lty="dashed", col="green")

# 6. Now for both NList4 and List8 combined.
# First get the mean durations for both lists in a convenient form:

Dmeans1 <- matrix(nrow=10,ncol=7)

for (pos in 1:7){
  for(digit in 0:9){
    Dmeans1[digit+1,pos] <- mean(c(Durs4[Nums4[,pos]==digit,pos], Durs8[Nums8[,pos]==digit,pos]))
  }
}

# ... and plot them:

cols <- c("red","darkgreen", "blue", "orange", "brown", "black", "maroon", 
           "purple", "rosybrown", "magenta" )
ltys <- c(1,2,3,4,5,6,1,2,3,4)
yr <- range(Dmeans1)

plot(1:7, Dmeans1[1,],
       type="n", ylim=yr,
       xlab="Digit Position", ylab="Mean Duration (seconds)",
       main="Mean Duration for Digits 0-9 (both lists)")

for (n in 0:9){
   points(1:7, Dmeans1[n+1,],
       type="b", col=cols[n+1], lty=ltys[n+1], pch=as.character(n))
}

#
# Now plot the relationship between duration in position 2 and position 3,
#     digit by digit:

lims <- range(c(Dmeans1[,2], Dmeans1[,3]))

plot(c(.3,.3),c(.3,.3), ylim=lims, xlim=lims, type="n",
     xlab="Position 2 Duration", ylab="Position 3 Duration",
     main="Internal Pre-boundary Lengthening, By Digit\n(both lists)")

text(Dmeans1[,2],Dmeans1[,3], labels=0:9, col=cols)
abline(a=0,b=1, col="black", lty="dotted")

# Fit a line just to monosyllabic digits c(1,2,3,4,5,6,8,9) 

monosyl <- c(1,2,3,4,5,6,8,9)+1  # add 1 since digit 0 is row 1
p11 <- glm(Dmeans1[monosyl,3] ~ Dmeans1[monosyl,2])

# And plot it:
abline(p11[1], lty="dashed", col="green")

# Ditto for position 6 and position 7:

lims <- range(c(Dmeans1[,6], Dmeans1[,7]))

plot(c(.3,.3),c(.3,.3), ylim=lims, xlim=lims, type="n",
     xlab="Position 6 Duration", ylab="Position 7 Duration",
     main="Final Pre-boundary Lengthening, By Digit\n(both lists)")

text(Dmeans1[,6],Dmeans1[,7], labels=0:9, col=cols)
abline(a=0,b=1, col="black", lty="dotted")

# Add line for c(1,2,3,4,9)
#     i.e. those without the complications of 
#      - final unreleased stop (8)
#      - two syllables (0,7)
#      - extra lengthening due to final fricative (5,6)

which <- c(1,2,3,4,9)+1
p22 <- glm(Dmeans1[which,7] ~ Dmeans1[which,6])
abline(p22[1], lty="dashed", col="green")