Dismiss Notice
Join Physics Forums Today!
The friendliest, high quality science and math community on the planet! Everyone who loves science is here!

MATLAB Non linear ODE system matlab fourier analysis problem

  1. Feb 22, 2012 #1
    Hello all.

    I am studying a system and want to investigate how the frequency of y(3) varies under different conditions. However, my the fft I perform on it tells me the frequency is zero, which must be incorrect. I have tried a stack of things but can't see what the problem is. I'm relatively new to the matlab environment.

    The main code:

    close all; clear all;
    Time control and initial conditions
    timeinmin=600; % Time in minutes
    tfin=timeinmin*60; % Time in seconds
    Initial Levels and system Saturation
    totIKK=0.08; % Total IKKn at start
    totNFkB=0.08; % Total NFkB cytoplasmic at start
    T=1; % Saturation of the system
    yint=[totNFkB,0,0,0,0,0,totIKK,0,0,0,0,totIKK,totNFkB,T]; % Initial conditions
    [t,y]=ode45(@nfkb,[0 tfin],yint); % calls Model
    Extract Data from the ODE
    fs = 1; % The sampling freqency in Hz
    time = 0:1/fs:tfin; % Define time of interesrt
    Datapoints = tfin*fs; % # of Data points being used
    yi = interp1(t,y(:,3),time);
    %plot(time/60,yi,'o',t/60,y(:,3)) % test to ensure it is working
    Fourier Analysis from
    NFFT = 2^nextpow2(tfin); % Next power of 2 from length of y
    Y = fft(yi,NFFT)/Datapoints;
    f = fs/2*linspace(0,1,NFFT/2+1);
    % Plot single-sided amplitude spectrum.
    title('Single-Sided Amplitude Spectrum of y(t)')
    xlabel('Frequency (Hz)')
    ylabel('Amplitude - |Y(f)|')
    % Get value of frequency
    [B,IX] = sort(2*abs(Y(1:NFFT/2+1)));
    BFloor=0.1; %BFloor is the minimum amplitude value (ignore small values)
    Amplitudes=B(B>=BFloor) %find all amplitudes above the BFloor
    Frequencies=f(IX(1+end-numel(Amplitudes):end)) %frequency of the peaks

    The ODEs:

    function dydt=nfkb(t,y)

    The parameters:

    kv=3.3; % C:N ratio
    kp=0.0006; % Conversion rate of IKKi
    ka=0.004; %activation rate of IKK
    ki=0.003; % Usage rate of IKK
    kda=0.0045; % Deg rate of A20 protein
    kbA20=0.0018; % A20 concentration
    kc2=0.074; % phosopho of IkB
    kc3=0.37; % phospho of IkB-NFkB
    kdegpin=0.1; % splitting of Ikb-NFkB due to phospho
    ka1=0.5; % nuclear IkB-NFkB formation
    kd1=0.0005; % dissisociation of nuclear IkB-NFkB
    kdegf=0.0005; % Deg rate of nuclear IkB
    kdegc=0.000022; % dissisociation of cyto IkB-NFkB
    ki1=0.0026; % nuclear entry of NF-kB
    ke1f=52*10^(-6); % nuclear exit of NF-kB
    ke1c=0.01; % Transport of IkB-NFkB out of nuc.
    ki2=0.00067; % nuclear entry of IkB
    ke2=3.35*10^(-4); % Nuclear exit of IkB
    h=2; % hill coefficient of IkB transcription
    k=0.065; % dissasociation constant for IkB tanscription with NFkB
    kitria=1.4*10^(-7); % trans rate of IkB
    ktria=0.5; % translation of IkB
    kdegtia=0.0003; % Deg rate of IkB RNA
    kitra=1.4*10^(-7); % Transcription rate of A20 RNA
    ktra=0.5; % transcription rate of A20
    kdegta=0.00048; % deg rate of A20 RNA

    Thank you in advance for any insights.
  2. jcsd
Share this great discussion with others via Reddit, Google+, Twitter, or Facebook

Can you offer guidance or do you also need help?
Draft saved Draft deleted