I have a function of type

$$C(t)=0.31 \left(e^{-0.267t}-e^{-1.72t}\right)$$

Where $t$ is in hours

The function represents the evolution in a window of 24 hours of a drug blood concentration, after a dose of a drug is given at $t = 0$ to a patient. I would like to use Fourier transform, to identify th best sampling step to reconstruct the concentration curve. Because it is a human whose blood we are sampling here I would like to reduce the sampling as much as possible, with a reconstruction error of around 20%.

I am using Matlab FFT function to compute the function's Fourier transform and plot the spectrum.

However I am not sure if my frequency are in Hertz or Millihertz. Because I am sampling in hours, I am expecting lower frequencies (millihertz), but I am not sure of the logic.

Here is my Matlab Code

clear all; clc;close all;
%% Time-domain signal
Ts=0.25; % time discretization step
t=0:Ts:24 % in hours
C=0.31* (exp(-0.267*(t))-exp(-1.72*(t)));


title('Time domain signal');
xlabel('Time in hours ')

%Frequency domain signal Via Matlab FFT function
title('Frequency domain signal via FFT of C(t)');
xlabel('Frequency in Hz or Millihertz?')

Provided you make your function causal, for instance with $C(t) =0$ for $t<0$ you get a function whose continuous Fourier transform is well defined (as the sum of two one-sided decaying exponentials). The first step is to compute it, for sanity check.

Then, $C$ is not band-limited, so a direct sampling (for Matlab) is likely to induce errors, but you can get a reasonable approximation, provided your sampling is fine enough. $0.25$ hour was a bit coarse, and on a smaller time range:

Time signal

I strongly recommand to tick your x-axis in seconds, to get sound frequencies in Herz. here I took:

Ts=0.1; % time discretization step
t=-0.03*3600:Ts:0.03*3600 % in hours (0.03) times seconds/hour
C=0.31* (exp(-0.267*(t))-exp(-1.72*(t)));
C(t<0) = 0;

Now, using a standard FFT code, you can get a 'reasonable' idea of the amplitude spectrum. Here I used FFTR.m, a code I made for a quick display of one half of the spectrum for a real function:


Time signal spectrum

What could be interesting for you is to superimpose the actual continuous Fourier transform, to check the differences.

  • $\begingroup$ Thank your for your answer, I know that at 0.25 hour/24hr i was forcing it a bit, the thing is the function that I am trying to sample represents the evolution of a drug blood concentration after a single dose is given to a patient & I would like to sample it in a way that I am able to reconstruct the curve and compute the area under the curve with limited aliasing. So may be I should have added these details into my question. Knowing that do you know of a way I can reduce my sampling window, with 4.5 Hz as the highest frequency per nyquist I am looking at sampling a human every 0.22 seconds $\endgroup$
    – Wazaa
    Jan 21 '16 at 19:18
  • $\begingroup$ Which is practically impossible in real life, I thought of decimation as a solution, but don't think it might work $\endgroup$
    – Wazaa
    Jan 21 '16 at 19:19
  • $\begingroup$ Knowing your constraints and a model of your data should certainly help, indeed $\endgroup$ Jan 21 '16 at 23:04
  • $\begingroup$ I have re-edited my question, please feel free to share any thoughts $\endgroup$
    – Wazaa
    Jan 22 '16 at 15:23

That function is not periodic. Using the Fourier transform to analyze its behaviour is questionable, to say the least. Maybe you're confusing things with the Laplace transform, which is more useful when analyzing control systems?

Why should the axis in your plot have Hertz as unit?

The discrete Fourier transform¹ does a base transform from a vector in $C^N$ to a vector in $C^N$; its physical significance is that the result is that each bin represents another frequency between $0$ and the the Nyquist rate (which, in the complex sampling case is the sample rate, and in the real sampling case is half the sample rate).

Since you only have a real function, I'd argue frequencies are relative to $\frac{f_{sample}}2$.

$f_{sample}$ is $\frac{4}{\text{hour}}=\frac4{60\cdot60\,\text{s}}= \frac1{900}\,\text{Hz}$ in your case. So each bin represents a frequency $$f\in [-\frac{f_{sample}}2 ; \frac{f_{sample}}2[\,=\,[-\frac{1}{1800}\,\text{Hz}; \frac{1}{1800}\,\text{Hz}[$$.

Now, frequency base vectors are uniformly spread over the nyquist rate, so the $n^\text{th}$ bin represents frequency $\frac{\frac{1}{900}\,\text{Hz}}{n_\text{bins}}$. $n_\text{bins}$ is the number of bins, which is the same as the input vector size; which is $4\cdot24=96$. So each bin is $\frac1{96\cdot900}\,\text{Hz}= 11.5741\,\mu\text{Hz}$. Not very handy.

I think you're somewhat confused what an FFT is. I think a textbook on signal theory will do you very well!

¹ you're using fft, which is a fast implementation of the DFT; not the continuous Fourier transform. That's an important difference.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.