# Filter design: MATLAB firls (Least-squares linear-phase FIR filter design) alternative?

**URL:** <https://discourse.julialang.org/t/filter-design-matlab-firls-least-squares-linear-phase-fir-filter-design-alternative/45249>\
**Category:** Signal and Image Processing\
**Tags:** dsp\
**Created:** [August 20, 2020, 7:44am UTC](https://discourse.julialang.org/t/filter-design-matlab-firls-least-squares-linear-phase-fir-filter-design-alternative/45249 "2020-08-20T07:44:16Z")\
**Posts on this page:** 18\
**Page:** 1

<div class="post-metadata">

**Author:** ![ElectronicTeaCup](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/electronicteacup/32/12647_2.png) [@ElectronicTeaCup](https://discourse.julialang.org/u/ElectronicTeaCup)\
**Post date:** [August 20, 2020, 7:44am UTC](https://discourse.julialang.org/t/filter-design-matlab-firls-least-squares-linear-phase-fir-filter-design-alternative/45249/1 "2020-08-20T07:44:17Z")

</div>

Hi! I am going through a MATLAB based [textbook](https://mitpress.mit.edu/books/analyzing-neural-time-series-data) for time frequency analysis of Neural data (MEG/EEG time series data) and got to a chapter which goes through filter construction. It’s primarily going over the usage of the functions `firls` (Least-squares linear-phase FIR filter design) and `fir1` (Window-based FIR filter design) (see: [docs](https://de.mathworks.com/help/signal/ref/firls.html) and more [docs](https://de.mathworks.com/help/signal/ref/fir1.html)). As a newcomer to this, I couldn’t find similar functions and was wondering if I could be pointed into the right direction. I have seen the filter docs in DSP.jl but they are a little difficult for me to understand. So

1. Which functions should I be using in the Julia ecosystem to do the same thing.
2. If the approach is different: how do they differ from Least-squares linear-phase FIR filter design and Window-based FIR filter design.

I should add I am looking for this because I’m doing the material in Julia, not in MATLAB. So I’m interested in learning and implementing, not comparing the methods between languages or anything.

---

<div class="post-metadata">

**Author:** ![ElectronicTeaCup](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/electronicteacup/32/12647_2.png) [@ElectronicTeaCup](https://discourse.julialang.org/u/ElectronicTeaCup)\
**Post date:** [August 20, 2020, 7:44am UTC](https://discourse.julialang.org/t/filter-design-matlab-firls-least-squares-linear-phase-fir-filter-design-alternative/45249/2 "2020-08-20T07:44:26Z")

</div>

(Sorry for the direct tag @ararslan, but I was informed you might be able to help with this, and I wanted to make a thread rather than a direct message.)

---

<div class="post-metadata">

**Author:** ![mbaz](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mbaz/32/17295_2.png) [@mbaz](https://discourse.julialang.org/u/mbaz)\
**Post date:** [August 20, 2020, 3:14pm UTC](https://discourse.julialang.org/t/filter-design-matlab-firls-least-squares-linear-phase-fir-filter-design-alternative/45249/3 "2020-08-20T15:14:27Z")

</div>

DSP.jl provides several ways to design filters. In your application, I suspect that the specific filtter design method is not very important. What you need to know is: the filter type (low pass, band pass, high pass); the sampling frequency; the filter order; and the cutoff frequency.

To design a filter, DSP.jl requires the specification of a design method, in addition to the parameters given above. For example, let’s say you want a low-pass filter with sampling frequency `200`, cutoff frequency `50`, and order `21`. For low pass filters, Butterworth is a good, easy design method. Then, you can run:

```julia
using DSP
fs = 200; # sampling frequency
fc = 50; # cutoff frequency
order = 21; # filter order
filter = digitalfilter(Lowpass(fc, fs=fs), Butterworth(order));

```

This produces a filter with the following magnitude response:  
 ![image](https://global.discourse-cdn.com/julialang/original/3X/9/1/9147c5c3857bb254815f61525f31adf121653f3e.png)

Now let’s say you want to filter a signal `x`. All you need to do is run

```
filtsignal = filt(filter, x);

```

---

<div class="post-metadata">

**Author:** ![ufechner7](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/ufechner7/32/51363_2.png) [@ufechner7](https://discourse.julialang.org/u/ufechner7)\
**Post date:** [August 22, 2020, 7:34am UTC](https://discourse.julialang.org/t/filter-design-matlab-firls-least-squares-linear-phase-fir-filter-design-alternative/45249/4 "2020-08-22T07:34:56Z")

</div>

But Butterworth is not linear phase. Try Bessel instead. See: [https://www.nuhertz.com/response/iir-and-analog-filters-basic-types/bessel-and-linear-phase-filters](https://www.nuhertz.com/response/iir-and-analog-filters-basic-types/bessel-and-linear-phase-filters)

---

<div class="post-metadata">

**Author:** ![mbaz](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mbaz/32/17295_2.png) [@mbaz](https://discourse.julialang.org/u/mbaz)\
**Post date:** [August 22, 2020, 2:22pm UTC](https://discourse.julialang.org/t/filter-design-matlab-firls-least-squares-linear-phase-fir-filter-design-alternative/45249/5 "2020-08-22T14:22:55Z")

</div>

You are correct, of course: indeed, the phase is slightly non-linear (the red line is the cutoff frequency):  
 ![image](https://global.discourse-cdn.com/julialang/original/3X/0/d/0d5df895a744e416eb1185cbb3c9a9a862c580ba.png)

Unfortunately DSP.jl does not support Bessel filters AFAICT. If perfectly linear phase is required, probably the best option is to use `remez`. Here’s a similar low-pass filter. Note that a larger order is required, and it has bandpass ripple (as often in filter design, it’s a matter of tradeoffs):

```julia
filter = remez(101, [(0, 50) => 1, (55, 100) => 0], Hz = 200);

```

![image](https://global.discourse-cdn.com/julialang/original/3X/8/9/8918944896d24d30797c66bcf3f8abc5e0907124.png)

---

<div class="post-metadata">

**Author:** ![ElectronicTeaCup](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/electronicteacup/32/12647_2.png) [@ElectronicTeaCup](https://discourse.julialang.org/u/ElectronicTeaCup)\
**Post date:** [August 26, 2020, 7:01pm UTC](https://discourse.julialang.org/t/filter-design-matlab-firls-least-squares-linear-phase-fir-filter-design-alternative/45249/6 "2020-08-26T19:01:19Z")

</div>

> [@mbaz](#):
>
> hat you need to know is: the filter type

Thank you for your replies! They really helped me getting started. One of the things that author does is define a specific shape of the filter (I think to make a smoother transition?) with the following Matlab code:

```MATLAB
nyquist = EEG.srate/2;
lower_filter_bound = 4; % Hz
upper_filter_bound = 10; % Hz
transition_width = 0.2;
filter_order = round(3*(EEG.srate/lower_filter_bound));

% create the filter shape

ffrequencies = [0 (1-transition_width)*lower_filter_bound lower_filter_bound upper_filter_bound (1+transition_width)*upper_filter_bound nyquist]/nyquist;
idealresponse = [0 0 1 1 0 0];
filterweights = firls(filter_order,ffrequencies,idealresponse);

```

So `ffrequencies` defines the frequencies of interest, normalized to the nyquist frequency. And then the shape of the filter is defined by `idealresponse`. Is there a way to make a smoother transition (assuming this is the reason) like this in Julia?

Also the filter order here comes out to something like 192, this doesn’t work with `Butterworth(192)`. What exactly is this value? And is it the same as the filter order for butterworth and the “order- `n` FIR filter” (sorry if this latter part may be beyond the scope of this discussion)

---

<div class="post-metadata">

**Author:** ![mbaz](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mbaz/32/17295_2.png) [@mbaz](https://discourse.julialang.org/u/mbaz)\
**Post date:** [August 26, 2020, 8:56pm UTC](https://discourse.julialang.org/t/filter-design-matlab-firls-least-squares-linear-phase-fir-filter-design-alternative/45249/7 "2020-08-26T20:56:59Z")

</div>

@ElectronicTeaCup That code can be adapted to work with DSP.jl’s `remez`, producing a very similar filter. See my example above, and follow up if you need more help.

---

<div class="post-metadata">

**Author:** ![ElectronicTeaCup](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/electronicteacup/32/12647_2.png) [@ElectronicTeaCup](https://discourse.julialang.org/u/ElectronicTeaCup)\
**Post date:** [September 2, 2020, 3:40pm UTC](https://discourse.julialang.org/t/filter-design-matlab-firls-least-squares-linear-phase-fir-filter-design-alternative/45249/8 "2020-09-02T15:40:16Z")

</div>

Thank you for the code and information. I checked out the documentation, and did the following to reproduce something similar to the MATLAB code I was using.

```julia
center_freq = 20 # in Hz
filter_frequency_spread = 6 # Hz +/- the center frequency
transition_width = 0.2
ffrequencies = [ 
    0, 
    (1-transition_width)*(center_freq-filter_frequency_spread), 
    (center_freq-filter_frequency_spread),
    (center_freq+filter_frequency_spread), 
    (1+transition_width)*(center_freq+filter_frequency_spread), 
    nyquist, 
] /nyquist
idealresponse = [0, 0, 1, 1, 0, 0,]

# Order hard set to 200 as in MATLAB code
filterweights = remez(
    200, 
    [
        (ffrequencies[1], ffrequencies[2]) => 0, 
        (ffrequencies[3], ffrequencies[4]) => 1, 
        (ffrequencies[5], ffrequencies[6]) => 0 
    ], 
    Hz = EEG["srate"]/nyquist
)

```

And then to take a look at what’s going on with a `fft`

```julia
filterweights = (firls(200,ffrequencies,idealresponse));

# compute its power spectrum
fft_filtkern = abs.(fft(filterweights))
fft_filtkern = fft_filtkern./maximum(fft_filtkern) # normalized to 1.0 for visual comparison ease

hz_filtkern = range(0,nyquist, length=101); # list of frequencies in Hz corresponding to filter kernel

```

 ![image](https://global.discourse-cdn.com/julialang/original/3X/b/d/bd45ed4470d8010c744a92d5fbb9666252c94e2e.png)

From the original MATLAB code:  
 ![image](https://global.discourse-cdn.com/julialang/original/3X/2/e/2ec544ddb3b47687a7003845e8c03fd226347714.png)

---

<div class="post-metadata">

**Author:** ![ElectronicTeaCup](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/electronicteacup/32/12647_2.png) [@ElectronicTeaCup](https://discourse.julialang.org/u/ElectronicTeaCup)\
**Post date:** [September 2, 2020, 6:57pm UTC](https://discourse.julialang.org/t/filter-design-matlab-firls-least-squares-linear-phase-fir-filter-design-alternative/45249/9 "2020-09-02T18:57:28Z")

</div>

When I change the `filter_frequency_spread` to 10, I get this (black is the `fft_filtkern` as in the example above:

 ![image](https://global.discourse-cdn.com/julialang/original/3X/9/3/9337d17df63445655854c436134a2da8e6b9356f.png)

I have explicity stated my desired response in `remez`, what am I doing incorrectly here?

---

<div class="post-metadata">

**Author:** ![mbaz](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mbaz/32/17295_2.png) [@mbaz](https://discourse.julialang.org/u/mbaz)\
**Post date:** [September 2, 2020, 9:58pm UTC](https://discourse.julialang.org/t/filter-design-matlab-firls-least-squares-linear-phase-fir-filter-design-alternative/45249/10 "2020-09-02T21:58:17Z")

</div>

The problem you’re seeing is that the filter you want cannot be created within the parameters you’ve specified. In broad terms:

- Narrow filters (relative to the Nyquist frequency), require large orders and/or wide transition bands.
- Very small transition bands also require large filter orders.
- Bandpass filters with lower cutoff frequency close to 0, or higher cutoff frequency close to Nyquist, are in general also hard to design.

The degrees of freedom you have available are: filter bandwidth, filter order, and transistion bandwidth. Relaxing one can improve the others. You may need to try a few combinations before you get a filter you like.

I gave it a shot with the bandwidth you need (20 Hz). I relaxed the transition bandwidth to 2 Hz and increased the order to 300. I also asked `remez` to work harder by increasing `maxiter`. This is what I got:

```julia
using DSP
fs = 250.0 # sampling frequency
fn = fs/2.0 # Nyquist frequency
B = 20.0 # filter bandwidth
fc = 20.0 # center frequency
fl = fc - B/2. # lower cutoff
fh = fc + B/2. # higher cutoff
tr = 2.0 # transition bandwidth
order = 400 # filter order

# calculate coefficients
coeff = remez(order, [(0, fl-tr) => 0,
                      (fl, fh) => 1,
                      (fh+tr, fn) => 0],
              Hz = fs,
              maxiter = 50)

# find and plot magnitude response
f = range(0, fn, length = 250)
fz = freqz(PolynomialRatio(coeff,[1.0]), f, fs)
mag = abs.(fz)
plot(f, mag)

```

![image](https://global.discourse-cdn.com/julialang/original/3X/d/0/d03014044529acd19998eece38b10976d067c0fc.png)

---

<div class="post-metadata">

**Author:** ![ElectronicTeaCup](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/electronicteacup/32/12647_2.png) [@ElectronicTeaCup](https://discourse.julialang.org/u/ElectronicTeaCup)\
**Post date:** [September 4, 2020, 1:51pm UTC](https://discourse.julialang.org/t/filter-design-matlab-firls-least-squares-linear-phase-fir-filter-design-alternative/45249/11 "2020-09-04T13:51:10Z")

</div>

Ah right, it’s majorly the asymetric transition zones in my example that seems to be causing the issue. Thank you so much, I didn’t have a clue that these parameters would make all the difference!

```julia
coeff = remez(order, [(0, 8) => 0,
                      (10, 30) => 1,
                      (36, fn) => 0],
              Hz = fs,
              maxiter = 50)

```

 ![image](https://global.discourse-cdn.com/julialang/original/3X/e/9/e98e94320d00dcd9f4a4260df583e029b9594378.png)

When I make the same transition on both sides (like your example) it works just fine

```julia
coeff = remez(order, [(0, 4) => 0,
                      (10, 30) => 1,
                      (36, fn) => 0],
              Hz = fs,
              maxiter = 50)

```

 ![image](https://global.discourse-cdn.com/julialang/original/3X/8/2/8269ffaf3e68cdbb94bac22c9327e7c9470535f9.png)

I do not know why exactly the transition zones differ on each end in the text I’m following, but I want to try to get the same outputs as the rest of the analysis. Do you think it would be possible to get a transition like the transition regions of the first block of code?

---

<div class="post-metadata">

**Author:** ![mbaz](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mbaz/32/17295_2.png) [@mbaz](https://discourse.julialang.org/u/mbaz)\
**Post date:** [September 4, 2020, 2:46pm UTC](https://discourse.julialang.org/t/filter-design-matlab-firls-least-squares-linear-phase-fir-filter-design-alternative/45249/12 "2020-09-04T14:46:00Z")

</div>

In my opinion, you’re hitting the limits of what `remez` can do. Unfortunately, `firls` and other more modern design algorithms don’t exist in DSP.jl yet, and `FIRWindow` does not allow specifying the transition bands.

Note that Matlab also struggles with this filter. Using `firls`:

```
h = firls(299, [0 8 10 30 36 fn]./fn, [0 0 1 1 0 0]);

```

 ![image](https://global.discourse-cdn.com/julialang/original/3X/f/3/f384ca889c4f0dc7c283d37215e61f8144617e09.png)

`fir1` does not support specifying transition bands (this command is equivalent to DSP.jl’s `FIRWindow` design method). The best results are obtained with `fir2`, but note that the transition bands seem quite a bit wider than specified:

```
h=fir2(299, [0 8 10 30 36 fn]./fn, [0 0 1 1 0 0]);

```

 ![image](https://global.discourse-cdn.com/julialang/original/3X/8/5/8589e71719a2035761a0e50d3eb2c2d7237b8eea.png)

What I would do in your case is to use `remez` to design a filter with narrow transition bands, flat passband, and good (40 dB or more) rejection, and move on 🙂 If you absolutely want to use the same filter as the textbook, you can generate the filter coefficients in Matlab and copy/paste them to your Julia program. If you don’t have access to Matlab, send me the code and I’ll send you the coefficients.

I hope I can find the time soon to implement `firls` and other design algorithms for DSP.jl, if nobody else beats me to it…

---

<div class="post-metadata">

**Author:** ![ufechner7](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/ufechner7/32/51363_2.png) [@ufechner7](https://discourse.julialang.org/u/ufechner7)\
**Post date:** [September 4, 2020, 3:00pm UTC](https://discourse.julialang.org/t/filter-design-matlab-firls-least-squares-linear-phase-fir-filter-design-alternative/45249/13 "2020-09-04T15:00:05Z")

</div>

> [@mbaz](#):
>
> `firls`

Scipy has it: [scipy.signal.firls — SciPy v1.11.4 Manual](https://docs.scipy.org/doc/scipy/reference/generated/scipy.signal.firls.html)

---

<div class="post-metadata">

**Author:** ![mbaz](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mbaz/32/17295_2.png) [@mbaz](https://discourse.julialang.org/u/mbaz)\
**Post date:** [September 4, 2020, 3:01pm UTC](https://discourse.julialang.org/t/filter-design-matlab-firls-least-squares-linear-phase-fir-filter-design-alternative/45249/14 "2020-09-04T15:01:32Z")

</div>

Right – it shouldn’t be hard to port it to Julia.

---

<div class="post-metadata">

**Author:** ![ElectronicTeaCup](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/electronicteacup/32/12647_2.png) [@ElectronicTeaCup](https://discourse.julialang.org/u/ElectronicTeaCup)\
**Post date:** [September 9, 2020, 4:55pm UTC](https://discourse.julialang.org/t/filter-design-matlab-firls-least-squares-linear-phase-fir-filter-design-alternative/45249/15 "2020-09-09T16:55:32Z")

</div>

I think I’ll follow your advice, and move on working with what I’ve got—thanks for the MATLAB offer btw. I really don’t want to go down the using MATLAB with Julia path—since the whole point was to learn Julia in this context. Thanks for all the support! I hope you end up porting `firls` 😀

---

<div class="post-metadata">

**Author:** ![mbaz](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mbaz/32/17295_2.png) [@mbaz](https://discourse.julialang.org/u/mbaz)\
**Post date:** [September 9, 2020, 4:56pm UTC](https://discourse.julialang.org/t/filter-design-matlab-firls-least-squares-linear-phase-fir-filter-design-alternative/45249/16 "2020-09-09T16:56:54Z")

</div>

You’re very welcome!

---

<div class="post-metadata">

**Author:** ![ElectronicTeaCup](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/electronicteacup/32/12647_2.png) [@ElectronicTeaCup](https://discourse.julialang.org/u/ElectronicTeaCup)\
**Post date:** [September 18, 2020, 2:49pm UTC](https://discourse.julialang.org/t/filter-design-matlab-firls-least-squares-linear-phase-fir-filter-design-alternative/45249/17 "2020-09-18T14:49:41Z")

</div>

I also just realized how easy it is to call this function using Pycall as an _alternative_ solution:

```julia
filty = pyimport("scipy")
firls = filty.signal.firls

```

Followed by code similar to the MATLAB approach:

```julia
center_freq = 20
filter_frequency_spread_wide = 10
ffrequencies = [ 
    0, 
    (1-transition_width)*(center_freq-filter_frequency_spread_wide), 
    (center_freq-filter_frequency_spread_wide), 
    (center_freq+filter_frequency_spread_wide), 
    (1+transition_width)*(center_freq+filter_frequency_spread_wide), 
    nyquist, 
]/nyquist
idealresponse = [0 0 1 1 0 0]
filterweightsW = (firls(201,ffrequencies,idealresponse))
plot((ffrequencies*nyquist),idealresponse[:])
fft_filtkern = abs.(fft(filterweightsW))
fft_filtkern = fft_filtkern./maximum(fft_filtkern) # normalized to 1.0 for visual comparison ease
hz_filtkern = range(0,nyquist, length=101)
plot!(hz_filtkern,fft_filtkern[1:ceil(Int,length(fft_filtkern)/2)],color=:black)

```

 ![image](https://global.discourse-cdn.com/julialang/original/3X/c/8/c8eee02aaf6fd1fdc127d8ac88cd00f526a20ea3.png)

---

<div class="post-metadata">

**Author:** ![LeoKoo](https://avatars.discourse-cdn.com/v4/letter/l/ccd318/32.png) [@LeoKoo](https://discourse.julialang.org/u/LeoKoo)\
**Post date:** [April 22, 2021, 6:49pm UTC](https://discourse.julialang.org/t/filter-design-matlab-firls-least-squares-linear-phase-fir-filter-design-alternative/45249/18 "2021-04-22T18:49:53Z")

</div>

Hi,

I packaged my own FIRLS implementation in Julia, see [this announcement](https://discourse.julialang.org/t/ann-firlsfilterdesign-jl/59820). Maybe you could give it a try and let me know what you think.

Cheers!