--- title: "Articulatory dynamics of (de)gemination in Dutch" author: "Patrycja Strycharczuk & Koen Sebregts" date: "20/02/2018" output: html_document --- ```{r setup, include=FALSE} knitr::opts_chunk$set(echo = TRUE) library(dplyr) library(lme4) library(effects) library(ggplot2) library(devtools) #install_github("kassambara/easyGgplot2") library(easyGgplot2) library(gamm4) library(mgcv) library(itsadug) #devtools::install_github("stefanocoretta/tidymv", build_vignettes = TRUE) library(tidymv) source_url("http://phon.chass.ncsu.edu/manual/tongue_ssanova.r") ``` # Duration analysis ```{r read-dur_data} load(file='data/dur.Rda') ``` ```{r dur-plot} dur.lmer=lmer(VR_dur~(1+Context|speaker)+(1|prompt)+(1|Block)+Context, data=dur) dur.lmer2=lmer(VR_dur~(1+Context|speaker)+(1|prompt)+(1|Block)+Context+Vowel, data=dur) anova(dur.lmer, dur.lmer2) summary(dur.lmer2) ef <- as.data.frame(effect("Context", dur.lmer2)) ef$Context = factor(ef$Context) ef$Context = factor(ef$Context, levels=c("Vr#C","Vr#V","VrV","V#rV", "Vr#r")) plot1 <- ggplot(ef, aes(x=Context, y=fit)) + geom_point(size=2) + geom_errorbar(aes(ymin=lower, ymax=upper, width=0.2))+ scale_x_discrete('')+ scale_y_continuous('Duration (ms)', limits=c(80,280))+ ggtitle('/ar/ duration')+ theme_bw()+ theme(text = element_text(size=12) ) cluster=subset(dur, Context=='Vr#r'|Context=='V#rV'|Context=='Vr#C') durC.lmer=lmer(cluster_dur~(1+Context|speaker)+(1|prompt)+Context, data=cluster) durC.lmer2=lmer(cluster_dur~(1+Context|speaker)+(1|prompt)+Context+Vowel, data=cluster) anova(durC.lmer, durC.lmer2) summary(durC.lmer2) ef2 <- as.data.frame(effect("Context", durC.lmer2)) ef2$Context = factor(ef2$Context) ef2$Context = factor(ef2$Context, levels=c("Vr#C", "V#rV","Vr#r")) plot2=ggplot(ef2, aes(x=Context, y=fit)) + geom_point(size=2) + geom_errorbar(aes(ymin=lower, ymax=upper, width=0.2))+ scale_x_discrete('')+ scale_y_continuous('Duration (ms)', limits=c(80,290))+ ggtitle('/ar(+C)/ duration')+ theme_bw()+ theme(text = element_text(size=12) ) cairo_pdf(file='/Users/mfuxjps2/Google Drive/Strycharczuk_Sebregts2015/paper/resubmission2/Figure3.pdf', width=5, height=3) ggplot2.multiplot(plot1, plot2, cols=2) dev.off() ``` #SS-ANOVA plots ```{r spline_data_import} load(file='data/spline_data.Rda') ``` ```{r ss-anova, echo=FALSE} df <- data.frame(X=as.numeric(character()), word=character(), ss.Fit=as.numeric(character()), ss.cart.SE=as.numeric(character()), ss.upper.CI.X=as.numeric(character()), ss.upper.CI.Y=as.numeric(character()), ss.lower.CI.X=as.numeric(character()), ss.lower.CI.Y=as.numeric(character()), Speaker=character() ) l = split(ssinput, ssinput[,1:2]) for (i in 1:length(l) ) { an=polar.ssanova((l[[i]][3:6]), crop=T)$pol.cart an$Speaker= l[[i]]$Speaker[1] an$time_point= l[[i]]$time_point[1] df=rbind(df,an) } ``` #SS-ANOVA plots ```{r ss plots, echo=FALSE} colnames(df)[2] <- 'Context' df=subset(df, Context!='VrV' & Context!='Vr#V') df$Context=factor(df$Context) df$time_point=factor(df$time_point) df$time_point=relevel(df$time_point, 'offset') df$time_point=relevel(df$time_point, 'r') df$time_point=relevel(df$time_point, 'onset') levels(df$time_point)<-c('/a/-onset', 'r', '/r/-offset') speaker=subset(df, Speaker=='DM2') speaker=subset(speaker, ss.Fit>-20) m <- matrix(c(1,1,1,2,3,4,5,5,5,6,7,8,9,9,9,10,11,12,13,13,13,14,15,16,17,17,17),nrow = 9,ncol = 3,byrow = TRUE) layout(mat = m,heights = c(1,2,1,2,1,2,1,2,2)) par(mar = c(0,0,1,0)) plot.new() text(0,0.8,"DF1",cex=1.5,font=2) #par(pty='s') for(n in 1:3){ speaker=subset(df, Speaker=='DF1') tp=speaker[speaker$time_point==levels(speaker$time_point)[n],] plot(c(), c(), ylim=range(speaker$ss.Fit), xlim=range(speaker$X), bty="n", xlab="", ylab="",asp=1, axes=F, main=levels(speaker$time_point)[n]) for(i in (1:length(levels(speaker$Context)))){ subdata=tp[tp$Context==levels(speaker$Context)[i],] lines(ss.Fit~X, data=subdata, col=i,lty=i) polygon(c(subdata$ss.upper.CI.X, rev(subdata$ss.lower.CI.X)), c(subdata$ss.upper.CI.Y, rev(subdata$ss.lower.CI.Y)), col=adjustcolor(i, alpha=(i-.5)/10) , border=F) } } plot.new() text(0,0.8,"DF2",cex=1.5,font=2) #par(pty='s') for(n in 1:3){ speaker=subset(df, Speaker=='DF2') tp=speaker[speaker$time_point==levels(speaker$time_point)[n],] plot(c(), c(), ylim=range(speaker$ss.Fit), xlim=range(speaker$X), bty="n", xlab="", ylab="",asp=1, axes=F, main=levels(speaker$time_point)[n]) for(i in (1:length(levels(speaker$Context)))){ subdata=tp[tp$Context==levels(speaker$Context)[i],] lines(ss.Fit~X, data=subdata, col=i,lty=i) polygon(c(subdata$ss.upper.CI.X, rev(subdata$ss.lower.CI.X)), c(subdata$ss.upper.CI.Y, rev(subdata$ss.lower.CI.Y)), col=adjustcolor(i, alpha=(i-.5)/10) , border=F) } } plot.new() text(0,0.8,"DF3",cex=1.5,font=2) #par(pty='s') for(n in 1:3){ speaker=subset(df, Speaker=='DF3') tp=speaker[speaker$time_point==levels(speaker$time_point)[n],] plot(c(), c(), ylim=range(speaker$ss.Fit), xlim=range(speaker$X), bty="n", xlab="", ylab="",asp=1, axes=F, main=levels(speaker$time_point)[n]) for(i in (1:length(levels(speaker$Context)))){ subdata=tp[tp$Context==levels(speaker$Context)[i],] lines(ss.Fit~X, data=subdata, col=i,lty=i) polygon(c(subdata$ss.upper.CI.X, rev(subdata$ss.lower.CI.X)), c(subdata$ss.upper.CI.Y, rev(subdata$ss.lower.CI.Y)), col=adjustcolor(i, alpha=(i-.5)/10) , border=F) } } plot.new() text(0,0.8,"DF4",cex=1.5,font=2) #par(pty='s') for(n in 1:3){ speaker=subset(df, Speaker=='DF4') tp=speaker[speaker$time_point==levels(speaker$time_point)[n],] plot(c(), c(), ylim=range(speaker$ss.Fit), xlim=range(speaker$X), bty="n", xlab="", ylab="",asp=1, axes=F, main=levels(speaker$time_point)[n]) for(i in (1:length(levels(speaker$Context)))){ subdata=tp[tp$Context==levels(speaker$Context)[i],] lines(ss.Fit~X, data=subdata, col=i,lty=i) polygon(c(subdata$ss.upper.CI.X, rev(subdata$ss.lower.CI.X)), c(subdata$ss.upper.CI.Y, rev(subdata$ss.lower.CI.Y)), col=adjustcolor(i, alpha=(i-.5)/10) , border=F) } } par(mar = c(1,1,1,1)) plot(0, 0, type = "n", ann = F, axes = F) legend(x = "top",inset = 0, legend=levels(tp$Context), col=c(1:3), lty=c(1:3), lw=2, cex=1, horiz=TRUE) m <- matrix(c(1,1,1,2,3,4,5,5,5,6,7,8,9,9,9,10,11,12,13,13,13,14,15,16,17,17,17),nrow = 9,ncol = 3,byrow = TRUE) layout(mat = m,heights = c(1,2,1,2,1,2,1,2,2)) par(mar = c(0,0,1,0)) plot.new() text(0,0.8,"DF5",cex=1.5,font=2) #par(pty='s') for(n in 1:3){ speaker=subset(df, Speaker=='DF5') tp=speaker[speaker$time_point==levels(speaker$time_point)[n],] plot(c(), c(), ylim=range(speaker$ss.Fit), xlim=range(speaker$X), bty="n", xlab="", ylab="",asp=1, axes=F, main=levels(speaker$time_point)[n]) for(i in (1:length(levels(speaker$Context)))){ subdata=tp[tp$Context==levels(speaker$Context)[i],] lines(ss.Fit~X, data=subdata, col=i,lty=i) polygon(c(subdata$ss.upper.CI.X, rev(subdata$ss.lower.CI.X)), c(subdata$ss.upper.CI.Y, rev(subdata$ss.lower.CI.Y)), col=adjustcolor(i, alpha=(i-.5)/10) , border=F) } } plot.new() text(0,0.8,"DM2",cex=1.5,font=2) #par(pty='s') for(n in 1:3){ speaker=subset(df, Speaker=='DM2') tp=speaker[speaker$time_point==levels(speaker$time_point)[n],] plot(c(), c(), ylim=range(speaker$ss.Fit), xlim=range(speaker$X), bty="n", xlab="", ylab="",asp=1, axes=F, main=levels(speaker$time_point)[n]) for(i in (1:length(levels(speaker$Context)))){ subdata=tp[tp$Context==levels(speaker$Context)[i],] lines(ss.Fit~X, data=subdata, col=i,lty=i) polygon(c(subdata$ss.upper.CI.X, rev(subdata$ss.lower.CI.X)), c(subdata$ss.upper.CI.Y, rev(subdata$ss.lower.CI.Y)), col=adjustcolor(i, alpha=(i-.5)/10) , border=F) } } plot.new() text(0,0.8,"DM4",cex=1.5,font=2) #par(pty='s') for(n in 1:3){ speaker=subset(df, Speaker=='DM4') tp=speaker[speaker$time_point==levels(speaker$time_point)[n],] plot(c(), c(), ylim=range(speaker$ss.Fit), xlim=range(speaker$X), bty="n", xlab="", ylab="",asp=1, axes=F, main=levels(speaker$time_point)[n]) for(i in (1:length(levels(speaker$Context)))){ subdata=tp[tp$Context==levels(speaker$Context)[i],] lines(ss.Fit~X, data=subdata, col=i,lty=i) polygon(c(subdata$ss.upper.CI.X, rev(subdata$ss.lower.CI.X)), c(subdata$ss.upper.CI.Y, rev(subdata$ss.lower.CI.Y)), col=adjustcolor(i, alpha=(i-.5)/10) , border=F) } } plot.new() text(0,0.8,"DM5",cex=1.5,font=2) #par(pty='s') for(n in 1:3){ speaker=subset(df, Speaker=='DM5') tp=speaker[speaker$time_point==levels(speaker$time_point)[n],] plot(c(), c(), ylim=range(speaker$ss.Fit), xlim=range(speaker$X), bty="n", xlab="", ylab="",asp=1, axes=F, main=levels(speaker$time_point)[n]) for(i in (1:length(levels(speaker$Context)))){ subdata=tp[tp$Context==levels(speaker$Context)[i],] lines(ss.Fit~X, data=subdata, col=i,lty=i) polygon(c(subdata$ss.upper.CI.X, rev(subdata$ss.lower.CI.X)), c(subdata$ss.upper.CI.Y, rev(subdata$ss.lower.CI.Y)), col=adjustcolor(i, alpha=(i-.5)/10) , border=F) } } par(mar = c(1,1,1,1)) plot(0, 0, type = "n", ann = F, axes = F) legend(x = "top",inset = 0, legend=levels(tp$Context), col=c(1:3), lty=c(1:3), lw=2, cex=1, horiz=TRUE) ``` # GAMM modelling ```{r read-total} load(file='data/art_dynamics.Rda') ``` ```{r gammFull} ad2 <- filter(art_dynamics %>% arrange(Date, time_norm) %>% start_event(column = "time_norm", event = c("Date"))) # mark beginning of events ad2$Vowel.ordered <- as.ordered(ad2$Vowel) contrasts(ad2$Vowel.ordered) <- 'contr.treatment' ad2$Context.ordered <- as.ordered(ad2$Context) contrasts(ad2$Context.ordered) <- 'contr.treatment' r.gam <- bam( LD1 ~ Context.ordered + s(time_norm) + s(time_norm, by = Context.ordered, bs = "cr") + s(time_norm, Speaker, bs="fs", k=10)+ s(time_norm, Prompt, bs="fs", k=10), data = ad2, method = "ML" ) r.gam2 <- bam( LD1 ~ Context.ordered + Vowel.ordered + s(time_norm) + s(time_norm, by = Vowel.ordered, bs = "cr") + s(time_norm, by = Context.ordered, bs = "cr") + s(time_norm, Speaker, bs="fs", xt="cr", m=1, k=10)+ s(time_norm, Prompt, bs="fs", xt="cr", m=1, k=10), data = ad2, method = "ML" ) compareML(r.gam, r.gam2) ad2$SpeakerContext <- interaction(ad2$Speaker, ad2$Context.ordered) r.gam3 <- bam( LD1 ~ Context.ordered + Vowel.ordered + s(time_norm) + s(time_norm, by = Vowel.ordered, bs = "cr") + s(time_norm, by = Context.ordered, bs = "cr") + s(time_norm, SpeakerContext, bs="fs", xt="cr", m=1, k=10)+ s(time_norm, Prompt, bs="fs", xt="cr", m=1, k=10), data = ad2, method = "ML" ) compareML(r.gam2, r.gam3) r1 <- start_value_rho(r.gam3) r.gam.AR <- bam( LD1 ~ Context.ordered + Vowel.ordered + s(time_norm) + s(time_norm, by = Vowel.ordered, bs = "cr") + s(time_norm, by = Context.ordered, bs = "cr") + s(time_norm, SpeakerContext, bs="fs", xt="cr", m=1, k=10)+ s(time_norm, Prompt, bs="fs", xt="cr", m=1, k=10), data = ad2, method = "ML", rho = r1, AR.start = ad2$start.event ) acf_plot(resid(r.gam.AR), split_by = list(ad2$Date)) ``` #Plot GAMM results ```{r GAMM plot} plot_smooth(r.gam.AR, view = "time_norm", rm.ranef=T, plot_all= "Context.ordered", rug = FALSE, col=c(4,3,1,2,5), main='', xlab='Time (normalized)' , legend_plot_all = list(x=2,y=0), ylab='LD1') text(0.0, 7, pos=4, labels='Vr#C', cex=0.8, xpd=TRUE) text(0.0, 4.5, pos=4, labels='Vr#V', cex=0.8, xpd=TRUE) text(0.0, 1.6, pos=4, labels='Vr#rV', cex=0.8 , xpd=TRUE) text(0.0, -1.2, pos=4, labels='VrV', cex=0.8, xpd=TRUE) text(0.0, -4.2, pos=4, labels='V#rV', cex=0.8, xpd=TRUE) fg=subset(ad2, Context=='Vr#r') fg$SpeakerContext=factor(fg$SpeakerContext) fg$Prompt=factor(fg$Prompt) sel=levels(fg$Prompt) speaker=levels(fg$SpeakerContext) speaker2=levels(ad2$Speaker) #par(fmrow=c(4,2)) cairo_pdf(file='Figure7.pdf', width=6, height=8) m <- matrix(c(1,2,3,4,5,6,7,8,9,9),nrow = 5,ncol = 2,byrow = TRUE) layout(mat = m,heights = c(2,2,2,2,1)) par(mar = c(4,4,4,4)) for(i in 1:8){ plot_smooth(r.gam.AR, view = "time_norm", cond=list(Context.ordered="Vr#r", Prompt=sel[1], SpeakerContext = speaker[i] ), rug = FALSE, main=speaker2[i], xlab='Time (normalized)', ylim=c(-6,7), ylab='LD1', cex.main=1.25, cex.lab=1.2, cex.axis=1.2) plot_smooth(r.gam.AR, view = "time_norm", cond=list(Context.ordered="Vr#r", SpeakerContext = speaker[i], Prompt=sel[2] ), add=T, lty=2, rug = FALSE, main='', xlab='', ylim=c(-6,7), ylab='') } par(mar = c(1,1,1,1)) plot(0, 0, type = "n", ann = F, axes = F) legend(x = "top",inset = 0, legend=list('paar reizen', 'paar raden'), lty=c(1:2), lw=1, cex=1, horiz=TRUE) dev.off() ```