Kubke Lab:Research/ABR/Notebook/2013/10/29

From OpenWetWare

< Kubke Lab:Research | ABR | Notebook | 2013 | 10(Difference between revisions)
Jump to: navigation, search
(Feedback from Andy)
Current revision (21:33, 28 October 2013) (view source)
(Feedback from Andy)
 

Current revision

Hearing development in barn owls Main project page
Previous entry      Next entry

Contents


General Entries

  • Played around with filter() and butter() in library.signal MF Kubke 20:11, 28 October 2013 (EDT)

Feedback from Andy

  • Filter: Higher order is sharper cutoff. A sharp butterworth would be a 5th or 6th order. I wouldn't use any higher.
  • Try an fft to see where the low pass cutoff should be. If the filter order is too high you get ringing. [ see http://t.co/QfrGaDjvj2 ]
  • entered wrong sampling rate - sampling rate is actually 50 kHz - need to recalculate the filters based on that. They are off by a factor of 10, so what says lowpassed at 100 is actually at 1k.

Personal Entries

Fabiana

  • File 2013-10-29-MFK.Rmd

Trying to apply a low pass filter to ABR data

Using package signal

install.packages("signal")
## Error: trying to use CRAN without setting a mirror
library(signal)
## Warning: package 'signal' was built under R version 3.0.2
## Loading required package: MASS
## 
## Attaching package: 'signal'
## 
## The following object is masked from 'package:stats':
## 
##     filter, poly

Read test file

dir()
##  [1] "2013-10-27-MFK.Rmd"  "2013-10-27.html"     "2013-10-27.R"       
##  [4] "2013-10-27.txt"      "2013-10-28-MFK.html" "2013-10-28-MFK.md"  
##  [7] "2013-10-28-MFK.Rmd"  "2013-10-29-MFK.Rmd"  "Copy189L0A.ABR"     
## [10] "Copy224l0a.txt"      "Copy233L0B.ABR.log"  "Copy233L0B.ABR.txt" 
## [13] "figure"              "knitr-testr.html"    "knitr-testr.md"     
## [16] "knitr-testr.Rmd"     "kntR2.html"          "kntR2.md"           
## [19] "kntR2.Rmd"           "my first knitr.html" "my first knitr.md"  
## [22] "my first knitr.Rmd"  "mydata_new"          "mydata_new.txt"     
## [25] "readrawabr.R"        "Sandbox1.R"          "Sandbox1.txt"       
## [28] "signal.pdf"          "texput.log"
mydata <- read.table("Copy233L0B.ABR.txt", header = T)
head(mydata)  #passed test
##   msec left right spont  bin
## 1 0.00 1000  1000  1000 1000
## 2 0.02 1000  1000  1000 1000
## 3 0.04 1000  1000  1000 1000
## 4 0.06 1000  1000  1000 1000
## 5 0.08 1000  1000  1000 1000
## 6 0.10 1000  1000  1000 1000
write.table(mydata, "mydata_new.txt")
mydata_new <- read.table("mydata_new", header = T)
head(mydata_new)  #passed test
##   msec left right spont  bin
## 1 0.00 1000  1000  1000 1000
## 2 0.02 1000  1000  1000 1000
## 3 0.04 1000  1000  1000 1000
## 4 0.06 1000  1000  1000 1000
## 5 0.08 1000  1000  1000 1000
## 6 0.10 1000  1000  1000 1000
mydata_new$bincomp[1:1000] <- (mydata_new$bin[1:1000] - (mydata_new$left[1:1000] + 
    mydata_new$right[1:1000]))
head(mydata_new)
##   msec left right spont  bin bincomp
## 1 0.00 1000  1000  1000 1000   -1000
## 2 0.02 1000  1000  1000 1000   -1000
## 3 0.04 1000  1000  1000 1000   -1000
## 4 0.06 1000  1000  1000 1000   -1000
## 5 0.08 1000  1000  1000 1000   -1000
## 6 0.10 1000  1000  1000 1000   -1000
write.table(mydata_new, "mydata_new.txt")  #passed test
head(mydata_new)
##   msec left right spont  bin bincomp
## 1 0.00 1000  1000  1000 1000   -1000
## 2 0.02 1000  1000  1000 1000   -1000
## 3 0.04 1000  1000  1000 1000   -1000
## 4 0.06 1000  1000  1000 1000   -1000
## 5 0.08 1000  1000  1000 1000   -1000
## 6 0.10 1000  1000  1000 1000   -1000

Test plotting for comparison


plot(mydata_new[c(1, 2)], type = "l", main = "original left")

plot of chunk unnamed-chunk-3

Now - try to apply a low pass filterusing butter() to determine the properties of the filter.

bf <- butter(3, 0.1)
test1 <- filter(bf, mydata_new$left)

par(mfrow = c(1, 2))
plot(mydata_new[c(1, 2)], type = "l", main = "original left")
plot(test1, main = "butter(3, 0.1)")

plot of chunk unnamed-chunk-4

Test type = “low” vs type = “high”

bf <- butter(3, 0.1, type = "low")
test1 <- filter(bf, mydata_new$left)

par(mfrow = c(1, 2))
plot(test1, main = "butter type low")

bf <- butter(3, 0.1, type = "high")
test1 <- filter(bf, mydata_new$left)
plot(test1, main = "butter type high")

plot of chunk unnamed-chunk-5

(guess default is low then)

Try changing filter properties - change filter order (first argument in function)

par(mfrow = c(2, 2))
plot(mydata_new[c(1, 2)], type = "l", main = "original left")

bf <- butter(1, 0.1, type = "low")
test1 <- filter(bf, mydata_new$left)
plot(test1, main = "butter(1, 0.1)")

bf <- butter(3, 0.1, type = "low")
test1 <- filter(bf, mydata_new$left)
plot(test1, main = "butter(3, 0.1)")

bf <- butter(15, 0.1, type = "low")
test1 <- filter(bf, mydata_new$left)
plot(test1, main = "butter(30, 0.1)")

plot of chunk unnamed-chunk-6

Not really sure what the order is doing - ????

Stick with filter order == 3, change crirical frequency W. W-> [0-1] where 1 is the Nyquist frequency (½ of the sampling rate )

head(mydata_new)
##   msec left right spont  bin bincomp
## 1 0.00 1000  1000  1000 1000   -1000
## 2 0.02 1000  1000  1000 1000   -1000
## 3 0.04 1000  1000  1000 1000   -1000
## 4 0.06 1000  1000  1000 1000   -1000
## 5 0.08 1000  1000  1000 1000   -1000
## 6 0.10 1000  1000  1000 1000   -1000

sampling rate => every 0.2 msec ()

Looks like sampling rate at 5 kHz => Nyquist = 2.5 KHz

par(mfrow = c(2, 2))
plot(mydata_new[c(1, 2)], type = "l", main = "original left")

bf <- butter(3, 0.01, type = "low")
test1 <- filter(bf, mydata_new$left)
plot(test1, main = "butter(3, 0.01 (25 Hz)")

bf <- butter(3, 0.1, type = "low")
test1 <- filter(bf, mydata_new$left)
plot(test1, main = "butter(3, 0.1 (250 Hz))")

bf <- butter(3, 1, type = "low")
test1 <- filter(bf, mydata_new$left)
plot(test1, main = "butter(30, 1 (2.5KHz))")

plot of chunk unnamed-chunk-8

Try the first derivative on the 250 LP


bf <- butter(3, 0.1, type = "low")
mydata_new$leftlow <- filter(bf, mydata_new$left)

par(mfrow = c(1, 3))
plot(mydata_new$msec, mydata_new$left, type = "l", main = "original left")
plot(mydata_new$msec, mydata_new$leftlow, type = "l", main = "lowpassed left")
mydata_new$diffleft2[2:1000] <- diff(mydata_new$leftlow)/diff(mydata_new$msec)
plot(mydata_new$msec, mydata_new$diffleft2, ylim = c(-3000, 3000), type = "l", 
    main = "first diff on lowpassed left")

plot of chunk unnamed-chunk-9

First deriv stil looks noisy - try different filters


bf <- butter(3, 0.04, type = "low")
mydata_new$leftlow <- filter(bf, mydata_new$left)

par(mfrow = c(1, 3))
plot(mydata_new$msec, mydata_new$left, type = "l", main = "original left")
plot(mydata_new$msec, mydata_new$leftlow, type = "l", main = "lowpassed  100 left")
mydata_new$diffleft2[2:1000] <- diff(mydata_new$leftlow)/diff(mydata_new$msec)
plot(mydata_new$msec, mydata_new$diffleft2, ylim = c(-3000, 3000), type = "l", 
    main = "first diff on lowpassed left")

plot of chunk unnamed-chunk-10


par(mfrow = c(1, 1))

plot(mydata_new$msec, mydata_new$diffleft2, ylim = c(-3000, 3000), type = "l", 
    main = "first diff on lowpassed left")

plot of chunk unnamed-chunk-10

Mmmmmmm….

Andy

  • Enter content here

Oris

  • Enter content here


Personal tools