Removing transients in highpass filtering with MATLAB

Question

I've made a simple first order IIR highpass filter with a zero at z = 1 and a pole at z = 0.9. Its frequency response looks like this:

Now, I filter a DC signal using this filter. Here's the MATLAB code I use to do it:

b = [1 -1]; % Zero at z = 1
a = [1 -0.9]; %Pole at z = 0.9

figure(1)
freqz(b, a)

t = 1:100;
x(1:length(t)) = 1; % Constant function

y = filter(b, a, x);

figure(2)
plot(t, x)
xlabel('Time');
ylabel('Input Signal');

figure(3)
plot(t, y)
xlabel('Time');
ylabel('Output Signal');

As my filter is highpass, I expect the DC to become zero, or atleast become severely attenuated. However, the output signal I get looks like this:

From my understanding, this exponential output is a transient produced because I haven't set the initial conditions correctly. Sure enough, setting x[-1] = 1 solves the problem. However, this works only for this particular input DC signal. For any general input signal, how do I set the initial conditions so that transients aren't produced?

Edit : I'm aware that the filtfilt() function does forward-backward filtering with transient minimization, but I really want to port the filter to an embedded platform, so I need to understand how transient removal works. Thanks in advance for the help!

Edit 2 : As suggested by Kuba Ober, I tried setting x[-1] as the value that it actually should have been. It works fine for a DC input, but here's what happened for a sinusoidal input:

clc; clear all;
p = 0.9

a = [1 -p]
b = [1 -1]


n = 1:100; % Samples
f = 0.2; % Frequency in Hz
Fs = 10; % Sampling rate in samples per second
t = n/Fs; % Time axis

x = sin(2*pi*f*t);

% Filter with the appropriate initial conditions
y = filter(b, a, x, filtic(b, a, [], [sin(2*pi*f*0)])); 

figure(1)
plot(t, x)
xlabel('Time');
ylabel('Input Signal');

figure(2)
plot(t, y)
xlabel('Time');
ylabel('Output Signal');

Here's the input signal :

And here's the output :

The first peak is visibly smaller than the second, which indicates some transients being present. I'm not entirely sure about this, but I think the reason it doesn't work is because just setting x[-1] is not enough, I also need to set y[-1]. The problem here is that there's no way to find out what y[-1] actually should be.

Edit 3 : Let me provide a little more info on the problem I'm working on. I'm trying to use filters to remove noise from ECG (Electrocardiogram) signals in an embedded platform. Here's a typical ECG signal, after filtering:

Here's what an ECG signal looks before filtering:

Note the DC offset in the signal before filtering. For filtering, I need a notch filter to remove high frequency power line noise and a highpass filter to remove the DC and the low frequency "drifting" of the signal.

The filters I use need to be linear phase, since the time domain morphology of an ECG signal is very important for diagnosis. However, my filter doesn't need to be real-time, as I'm doing the processing offline after acquiring the ECG signal from the patient. So, for implementing nonlinear phase IIR filters, I'm currently using forward-backward zero phase filtering.

One opinion that's shared by @Matt L. and @Royi is that transients are unavoidable in real-time filtering and that I should use a longer input signal and crop off the first few seconds of the filtered output instead. This is something I'd like to avoid, as acquiring long ECGs from a living patient is somewhat difficult. Also, I do not have to filter in real-time, so any technique of transient removal that hinges on knowing the entire signal in advance is perfectly all right. Any help is appreciated!

Answer 1

Your first order filter recursion for some real constants $a,b,c$ is

$$ y[n] = a x[n] + b x[n-1] - c y[n-1] $$

with the two initial memory states $x[-1]$ and $y[-1]$ at $n=0$.

Your "no transient" condition can be translated to $y[0]=0$ and a necessary second condition so that you can solve for both of your memory states. The second condition could be, that the discrete derivative of $y$ also vanishes at $n=0$, so $y[0]-y[-1]=0$. You can also take any other condition that seems sensible to you.

The two equations give you a unique solution for the two unknown memories, namely: $$y[-1] = 0$$ and $$x[-1] = -\frac{a}{b}x[0] $$

Alternatively, your conditions may better be chosen as $y[0]=0$ and $y[0]-y[-1]=x[0]-x[-1]$ in order to capture the initial slope of the input. The resulting recursion equation at $n=0$ is then $$0=a x[0]+b x[-1]+c(x[0]-x[-1])$$ giving you the solution $$x[-1]=-\frac{a+c}{b-c}x[0]$$ and $$ y[-1]= -\left(1+\frac{a+c}{b-c}\right)x[0]$$

(Please check my calculations!)

But in general you cannot expect a simple initial condition to give you the same result as knowing the signal history. So you can only take this to a certain point and in general it would probably be better if you just discarded the transient response of your output.